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ABSTRACT 


'«r 


The  utility  of  a  computationally  simple  yet 
cry'ptologically  robust  rule  for  generating  pseudorandom 
bitstreams  cannot  be  overstated.  In  most  applications  we 
strive  to  detect  and  avoid  chaotic  behavior;  here  we  embrace 
a  particular  chaotic  discrete  dynamical  system  to  exploit  its 
use  as  a  driver  for  a  pseudorandom  number  generator.  The  map 
from  the  Henon  attractor  to  the  binary  domain  {0,1}  proposed 
by  Forre/Heyman  has  been  tested  cryptologically  and 
statistically  with  mixed  results.  In  this  thesis  we 
mathematically  evaluate  this  symbolic  dynamics  scheme  and 
investigate  more  rigorously  its  utility  as  a  pseudorandom 
number  generator.  Specifically,  we  demonstrate  how  the 
property  of  being  one-to-one  holds,  but  that  the  property  of 
being  onto  does  not . 
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I .  INTRODUCTION 


A.  HtNON'S  DISCRXTS  DYNAMICAL  SYSTSM 

Michel  Henon  (1931-  )  was  a  French  mathematician  and 

astronomer  who  was  interested  in  modeling  astronom.icai 
behavior.  Henon  learned  how  the  m.eteorologist -mathemat  ician 
Edward  N.  Lorenz  (1917-  )  had  developed  a  three-dimensional 

attractor  to  model  the  complicated  dynamics  of  thermal 
convection.  An  attractor  is  characterized  by  the  following 
definition . 

DEFINITION:  Let  V  be  a  subset  of  R*  and  F;  V  — >  R',  where 
m  =  1,2,3.  Let  A  be  a  subset  of  V.  Then  A  is  an  attractor  of 
F  if: 

i.  A  is  a  closed  invariant  subset  of  V  under  F.  The  set 

A  in  R'  is  invariant  under  F:  V  — >  R""  if  F(A)  c  A. 

li.  There  is  a  neighborhood  U  of  A  such  that  if  v  is  in  U, 

F'lv)  — ^  A  as  n  — »  oo  ,  where  F'  denotes  the  n-fold  composition  of 

F  [Ref .l:p. 147, 200]  . 

The  Lorenz  attractor  is  generated  by  a  system  of  three 
differential  equations  and  relies  upon  integration  which  does 
not  lend  itself  to  timely  and  accurate  computer  calculations. 
[Ref.  2]  Henon  sought  to  build  a  simple  discrete  dynamical 
system  that  would  retain  the  qualitative  properties  of  the 


continuous  Lorenz  system.  His  goal  was  achieved  through  study 
of  the  Poincare  section  of  the  Lorenz  equations.  [Ref. 3] 
Thus,  the  Henon  attractor  is,  to  some  extent,  a  two- 
dimensional  version  of  the  Lorenz  attractor.  The  Henon  map  is 
the  following  quadratic  recurrence  which  maps  the  Euclidean 
plane  into  itself.  [Ref. 7] 


H  :  R- 


X.. ;  =  1  -  axn  +  y. 
y,.:  =  bx. 


Each  point  in  R-  has  a  unique  image  and  preimage  under  H.  The 
samiS  is  true  of  the  inverse  of  the  Henon  map.  The  constants 
a  and  b  are  real  parameters.  We  can  write  the  map  as  a 
transformation  in  the  following  form; 

H(x,y)  =  (1  -  ax,-  +  y, ,  bx, )  . 

Understanding  the  dynamics  of  the  map  is  easier  if  we 
decompose  it  into  the  following  independent  transformations: 


H;(x,y)  =  (x,y  +  1  -  ax- 

H_(x,y)  =  (bx,y) 

H-.(x,y)  =  (y,x) 


such  that 


H  {x,y)  =  H,  (H;  (H;  (x,y) 


•  • 


long  as  Che  values  tend  to  the  attractor  and  the  first,  sa\-, 
ICC  Iterates  are  ignored.  The  first  points  generated  by  the 
T.ap  10, ay  not  be  close  to  the  attractor,  but  the  generated 
ctints  get  closer  and  closer  to  the  invariant  set. 

ZEF tl.'ITI CM :  Let  x  belong  to  the  domain  of  a  function  f;  x  is 

:f  period  n  if  f  ;x  =  x,  and  is  of  prime  period  n  if 

X,  f  ,.x  '  ,  f- ix;  ,  .  ■  .  .  f’  (x)  are  distinct.  If  x  has  perici  n  then 
the  orbit  {x, f (x) , f- (x)  ,  .  .  .  ,  f  ■  ;x) }  is  a  periodic  orbit 
;  .Fef  .  1 :  p  .  2 1  ]  .  If  n  -  1,  then  x  is  said  to  be  a  fixed  point. 


The  merit  of  Henon's  choice  of  ;a,b)  parameters  can  be 
found  by  observation  of  the  respective  bifurcation  diagrams. 

bifurcation  diagram  gives  us  information  regarding  the 
periodic  and  nonperiodic  orbits  cf  a  map  as  a  par-meter  is 
altered.  In  a  bifurcation  diagram,  the  verticil  ax-s 
represents  iterates  cf  the  variable  x  as  a  param.eter  is  varied 
on  the  horizontal  axis.  [Ref. rip. 52]  If  we  take  a  ver'-cal 
slice  at  a  particular  param.eter  we  get  an  laei  if  -he 
attractor's  structure  for  that  parameter.  In  Figure  2  the  h- 
parameter  is  held  fixed  at  0.3  and  the  a-paramieter  is  ••iried 
f  ror.  C  .  2  to  1.42. 
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Figure  2  Bifurcacion  Diagram;  b-paraunecer  fixed 


Notice  that  from  about  a  =  0.1  to  a  =  0.4  we  have  a  pome 
attractor,  that  is,  the  entire  attractor  consists  of  a  single 
fixed  point.  From  about  a  =0.4  to  a  =0.9  we  have  a  period- 
two  attractor;  that  is,  the  attractor  is  made  up  of  two 
points.  These  period-doubling  bifurcations  persist  as  the 
parameter  a  is  increased,  until  a  critical  value  of 
approximately  a  =  1.1  is  reached,  where  we  enter  the  chaotic 
regime.  Notice  that  Henon's  classical  a-parameter  1.4  lies  in 
such  an  area  of  wildly  complicated  dynamics,  sometimes  called 
a  chaotic  band  [Ref. 6]. 

In  Figure  3  the  a-parameter  is  held  fixed  at  1.4  and  the 
b-parameter  is  varied  from  0.0  to  0.32.  Again,  Henon's  b- 


< 
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value  0.3  not  surprisingly  lies  in  an  area  of  chaotic 
activity. 


1.3  r 


5  ■ _ I _ ■  _ _  _ ! 

0  0.05  0.1  0.15  0.2  0.25  0.3  0.35 

b-p«cameter 

Figure  3  Bifurcation  Diagraun:  a-parameter  fixed 

B.  NATURE  AND  STRUCTURE  OF  THE  H^ON  ATTRACTOR 

Because  the  Henon  transformation  is  quadratic,  it  is 
possible  that  even  relatively  small  (x-y,)  values  can  produce 
an  Iterated  sequence,  or  orbit,  that  escapes  to  infinity. 
[Ref .4 :p. 664]  The  Henon  attractor  has  an  associated  basin  of 
attraction  of  points  (found  numerically)  which  in  forward 
iteration  are  drawn  to  the  attractor.  A  point  in  the  basin  of 
attraction  thus  tends  to  the  invariant  set.  The  Henon 
attractor  also  has  an  associated  area  known  as  the  crapping 


•  • 
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region.  The  trapping  region  is  a  quadrilateral  within  the 
basin  of  attraction  defined  by  the  vertices: 

A  =  {-1.33,0.42)  B  =  (1.32,0.133) 

C  =  (1.245,-0.14)  D  =  (-1.06, -0.5) 

If  ;x-,y  )  is  chosen  from  this  region,  all  subsequent  iterates 
will  remain  in  the  region  [Ref . 5 : p . 664 ]  ;  hence  the  trapping 
region  is  invariant.  Let 

=  Henon  attractor 
=  Basin  of  attraction 
TR  =  Trapping  region. 

Then,  we  have  that 

H;,  C  TR  c  B,. 

This  relationship  is  shown  explicitly  (excluding  the  trapping 
region)  in  Figure  4;  the  shaded  region  is  the  basin  cf 
attraction . 
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Flgur*  4  Basin  of  attraction  (shaded 
region)  and  Henon  attractor.  This 
figure  is  from  reference  5,  p.665. 

Therefore, 

i.  if  (Xg.y^)  e  TR,  then  (x„,y„)  e  TR  for  all  n  >  0. 
ii.  if  e  but  (x^.y^)  e  TR,  then  (x.,,yj  may  not 

belong  to  TR  for  n  <  N.  However,  from  the  description 
of  the  basin  of  attraction  we  know  that  for  all  n 
sufficiently  large,  (x„,y„)  e  TR. 

The  b-parameter,  0.3,  is  in  the  diaslpatlve  term  of  the 
relation.  The  Jacobian  of  the  Henon  transformation  is  (-b), 
and  hence  the  map  is  dissipative  when  |b|  <  1.  [Ref . 1 :p . 168]  . 
Therefore,  should  we  iterate  the  entire  trapping  region  under 


the  transformation  for  b  =  0.3,  this  area  would  shrink  and  in 
the  limit  tend  to  zero.  The  Henon  attractor  thus  has  zero 
area.  Shown  symbolically, 

H.  =  liiTt  H'  (TR) 
or 

=  rv,  "  H"  (TR)  . 


In  fact,  any  point  belonging  to  the  trapping  region  can  be 
iterated  to  construct  a  dense  orbit  on  H^,  except  the  periodic 
points  which  are  countable  and  hence  have  Lebesgue  measure 
zero.  The  visual  result,  in  the  limit,  is  the  same  for  any 
initial  condition  (less  a  finite  number  of  initial  iterations 
to  assure  close  proximity  to  the  attractor) .  Two  different 
points  will  generate  different  sequences  but  the  limit  set, 
under  the  above  qualification,  will  appear  to  be  the  sa.me 
visually.  Furthermore,  in  general,  there  is  no  correlation 
between  the  trajectories  of  two  different  initial  points, 
except  possibly  in  the  first  relatively  few  iterations  should 
the  initial  points  be  neighbors.  (We  disregard  here  the 
instance  where,  say,  trajectory  1  has  common  points  with 
trajectory  2  but  for  different  iterates.  In  this  case  the 
orbits  are  shifted  versions  of  one  another.)  This  laclc  of 
correlation  is  related  to  sensitive  dependence  on  initial 
conditions . 
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The  Henon  equations  (1.1)  and  (1.2)  (with  the  classical 
parameters  (a,b)=(1.4.0.3) )  may  be  written  in  one  variable  as : 

x„,,  =  1  -  1.4x„^  f  0.3x„.,.  (1.3) 

There  are  two  fixed  point  solutions  to  the  recurrence,  where 

Xn*l  =  X„  =  X„.,. 

This  occurs  when  (x,y)  =  {. 6314 1894 ... )  (FP2)  and  (x,y) 

=  (-1.1314...,  -.3393...)  (FPl).  The  fixed  points,  trapping 
region  and  H4non  attractor  appear  in  Figure  5. 


Figure  5  Fixed  points  of  Henon  attractor,  trapping 
region  and  eigenvectors  associated  with  FP2 


The  left  fixed  point  FPl  is  repelling  and  is  noticeably  not  on 
the  attractor  or  even  within  the  trapping  region.  The 
remaining  fixed  point  FP2  is  a  saddle  point  and  belongs  to  the 
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attractor.  The  Jacobian  matrix  at  FP2  has  associated 
eigenvalues  X-i  and  X,,  with  slopes  Pi  and  p,  for  the  respective 
eigenvectors : 

X,  =  0.1559.  .  .  p..  =  1.9237.  .  . 

Xj  =  -1.9237. . .  p,  =  -0.1559.  .  . 

Figure  5  also  shows  us  the  strong  local  effect  of  the 
eigenvalues  at  FP2 .  Longitudinal  spreading  is  due  to  Xj . 
Since  [XjI  is  greater  than  1,  points  are  forced  outward  from 
FP2  along  the  longitudinal  axis  forming  a  roughly  broken  line, 
at  least  locally.  The  eigenvalue  Xi  causes  a  transverse 
contraction  of  points.  Figure  6  shows  a  detail  of  this 
structure.'  The  dense  transverse  structure  is  revealed  through 
local  magnification  in  Figure  7. 


Figure  6  Magnification  of  area  around  fixed  point 
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Figure  7  Decail  of  dense  transverse  structure  around 
fixed  point 
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As  magnification  is  increased  and  more  points  are  plotted  the 
attractor  begins  to  take  on  the  appearance  of  distinct  bands 
with  a  cross-section  that  has  the  nature  of  a  one-dimensional 
Cantor  set  --  a  totally  disconnected  set  of  points  with 
patterns  that  persist  on  infinitely  small  scales.  This  is  the 
self -similar  nature  of  the  attractor.  Thus,  at  least  on  a 
level  local  to  FP2  we  observe  a  fractal  quality  that  is  common 
in  chaotic  attractors  [Ref . 5 :p . 670]  . 

C .  CHAOS 

Because  the  field  of  chaos  is  in  its  relative  infancy  not 
all  terms  are  universal.  The  most  mathematically  rigorous  and 
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widely  accepted  definition  of  chaos  for  discrete  dynamical 
systems  is  the  following: 

DEFINITION:  Let  S  be  a  set.  A  mapping  f :  S  — ^  S  is  said  to  be 

chaotic  on  S  if 

1.  f  has  sensitive  dependence  on  initial  conditions, 

2.  f  is  topologically  transitive, 

3.  the  set  of  periodic  points  is  dense  in  S.  [Ref.8:p.50] 


* 


A’e  will  now  elaborate  on  each  of  these  concepts  in  a  bit  more 
detail  starting  with  sensitive  dependence  on  initial 
conditions . 


DEFINITION:  A  mapping  f:  S  -*  S  has  sensitive  dependence  on 
initial  conditions  if  at  every  x  €  S,  there  exists  an  e  >  0 
such  that  for  each  S  >  0,  there  is  a  y  in  S  and  a  positive 
integer  n  such  that 

|x-y  |<  5  and  IfCx)  -  f'{y)|>  £■  [Ref.l:p.84] 

Based  on  empirical  evidence,  the  Henon  attractor  is  widely 
believed  to  possess  sensitive  dependence  on  initial 
conditions.  [Ref. 7]  The  stretch-and-foid  phenomenon  mentioned 
previously  ensures  that  close  points  will  not  be  neighbors 
after  a  finite  number  of  iterations.  Figure  8  shows  a  plot  of 


•  • 
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the  Euclidean  distance  between  the  orbits  of  two  initial 


conditions . 


Figura  8  Sensitive  dependence  on  initial  conditions 


The  two  points  are  the  origin  (0,0)  and  a  point  (0,e')  where  ^ 

e'  is  just  within  machine  epsilon,  that  is  e'  <  e.  Machine 
epsilon  in  this  case  is  approximately  2.22  x  lO'-'’.  Within 
seventy  iterations  the  orbits  are  completely  different.  0 

Of  course,  since  both  initial  conditions  belong  to  the 
trapping  region  the  maximum  difference  between  the  orbits  is 
bounded  by  the  extremes  of  the  trapping  region.  Although  this  • 

is  only  the  first  of  three  criteria  for  chaos  some  accept  the 
dynamics  on  the  Henon  attractor  as  chaotic  based  on  sensitive 
dependence  alone.  [Ref . 4 :p . 671]  Sensitive  dependence  brings  • 
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out  an  interesting  paradox.  Globally  we  have  an  attractor 
with  respect  to  the  basin  of  attraction;  however,  locally  for 
the  points  of  the  invariant  set,  the  orbit  jumps  around  the 
attractor  chaotically,  seemingly  repelled  by  each  point . 

Topological  transitivity,  also  called  mixi.ng,  ensures  chat 
the  sec  on  which  the  dynamical  system  is  defined  cannot  be 
divided  into  two  disjoint  open  sets  which  do  not  “interact" 
ior  mere  precisely,  are  invariant)  under  the  system, 

DEFINITION:  A  mapping  f:  S -4 S  is  topologically  tranaitive  if 
for  any  pair  of  open  sets  U,  V,  which  are  subsets  of  S,  there 
exists  k  >  C  such  that  f'MU)  n  V  *  0.  [Ref  .8:p.45] 

Essentially,  for  a  system  to  be  topologically  transitive 
It  must  be  true  that  for  any  two  arbitrarily  small  regions  U 
and  V  (where  U  and  V  are  not  points)  we  can  find  an  initial 
point  X-  in  U  whose  orbit  will  enter  the  other  region  at  some 
iteration.  [Ref . 5 :p . 554] 

For  a  system  to  possess  dense  periodic  points  means  that 
in  an  arbitrarily  small  neighborhood  of  any  point  belonging  to 
the  attractor  there  must  be  a  periodic  point.  This  condition 
suggests  paradoxically  that  one  of  the  qualities  for  "chaos" 
IS  a  high  degree  of  "order". 

Although  there  exists  vast  numerical  evidence  to  support 
the  claim  that  Henon's  dynamical  system  possesses  all  three  of 
these  criteria  for  chaos,  none  of  the  criteria  have  been 
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rigorously  proven.  [Ref . 5 : p . 671 ]  [Ref . 19 :p . 152 ]  To  keep  this 
in  perspective,  it  must  be  acknowledged  that  our  definition  is 
more  stringent  than  proving  the  existence  of  a  dense  orbit  and 
a  positive  Lyapunov  exponent,  a  commonly  accepted  alternative 
definition.  [Ref. 14] 

D.  THE  SHADOWING  LEMMA  -  TRUSTING  THE  C(»(PUTER 

The  study  of  theoretical  chaotic  dynamical  systems  would 
be  extremely  handicapped  if  we  could  not  trust  our  computer 
representations  of  these  systems.  Because  a  computer  is 
limited  in  terms  of  its  accuracy  by  rounding  error,  it  is 
appropriate  for  us  to  ask  whether  our  computer  calculation  is 
an  approximation  of  some  true  orbit  of  the  system,.  [Ref.  9] 
The  Shadowing  Lemma  shows  that  we  can  trust  the  computer 
imiage . 


DEFINITION  1;  The  map  f :  R"  — »  R  is  a  class  C'  map  if  it  is  p 
times  continuously  differentiable  with  respect  to  x  ,  x_ ,  .  . 
for  some  p,  1  <  p  <  oo.  a  mapping  f  :R"— >  R’"  is  class  C"  if  each 
component  f  of  f, (where  i=l,...,m)  is  p  times  continuously 
differentiable  [Ref . 10 :p  .  1-2 ]  . 

DEFINITION  2:  A  mapping  f:R"->R*  is  a  dif f eoBiorphism  if  it  is 
one-to-one,  onto  '  and  both  f  and  f  ■  are  continuously 
differentiable  maps;  f  is  a  C-dif f eomorphism  if  both  f  and  f 
are  C"  maps  [Ref . 10 :p . 2 ]  . 
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DEFINITION  3:  Let  f:  R"  — »  R‘  be  a  dif f eomcrphism.  A  set  A  is 
said  to  be  a  hyperbolic  set  under  f  if: 

1.  for  all  points  p  belonging  to  A,  there  is  a  set  of 
lines  EMpi  and  EMp)  in  the  tangent  plane  at  p  which  are 
preserved  by  the  Jacobian  of  f  at  p.  (E-(p)  is  called  the 
stable  line  and  E(p)  is  called  the  unstable  line) . 

2.  E'-'p)  and  E'{p)  vary  continuously  with  p. 

1.  There  is  a  constant  X>1  such  that  |J(p) iv) |>X|v|  for 
all  veE  ( p  )  and  |  J  Mp )  ( v)  |  >X|  v|  for  all  veE®  (p )  .  [Ref  .  S  :  p  .  2  3  6  - 
22T] 

Suppose  f  satisfies  the  three  previous  definitions  for  some  p 
with  an  associated  hyperbolic  invariant  set  A.  The  true 
(forward)  orbit  of  a  point  z.  under  the  map  f  is  given  by  the 
sequence  {2;},-“  such  that  z,  =  f-(z-).  (Ref  .  10  :p .  351  ]  Due  to 
finite  precision  the  exact  calculation  of  the  sequence 
{flz.)}”  is  impossible.  Because  we  are  unable  to  obey  the 
function  z.  =  f-(z-)  in  order  to  find  an  exact  orbit,  let  us 
instead  define  a  pseudo-orbit  or  more  specifically  e-pseudo- 
orbit  . 


•  • 


DEFINITION  4:  An  e -pssudo-orbit  is  a  sequence  of  points, 
{y.}  ”  ,  such  that  y,  €  A,  and 

ci(y..,,  f  (y.)  )  <  e 
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where  d  is  a  metric  on  R".  [Ref  .  10 -.p  . 352  ]  The  e-pseudo-orbit 
is  pictured  at  the  top  of  Figure  9. 


Notice  that  we  adjust  each  iterate  by  dropping  f(y,)  and  then 
choosing  y^.l  to  be  within  e  of  f  (yj  .  We  can  see  that  the 
term  pseudo-orbit,  therefore,  is  appropriate.  One  final 
definition  is  necessary  for  us  to  understand  the  Shadowing 


Lemma . 


DEFINITION  5:  Suppose  {x_}-"  is  an  actual  orbit,  that  is,  x  eA 
and  X  =  f'(x  )  ;  {x  is  a  5-ahadow  of  the  pseudo-crbit  '.y  }  " 

d  'x  , y  )  <  5 . 

The  5-shadow  of  the  pseudo-orbit  is  shown  at  the  boctorr.  of 
Figure  9.  The  Shadowing  Lemma  assures  us  that  within  an 
arbitrarily  small  distance,  5,  of  an  £ -pseudo-orbit  there  is 
an  exact  orbit.  Stated  succinctly; 

Shadowing  Lemma:  If  Definition  2  holds  for  the  map  f  and  A  is 
a  hyperbolic  invariant  set  on  R",  then  for  every  5  >  C  there 
13  an  e  >  0  such  that  every  e -pseudo-orbit  in  A  is  8- 
shadowed  by  the  actual  orbit  of  some  point  x  e  A. 
[Ref . 10 :p. 352] 

Although  the  Shadowing  Lemma  is  very  powerful,  it  does  not 
explicitly  speak  to  us  about  the  computed  orbit,  ca;.l  it  w  , 
for  which  we  hope  to  find  to  an  exaot  orbit  nearby. 

Our  computed  orbit,  w  ,  is  not  the  true  orbit  •:  f  'y  '  }  “  we 
are  lucking  for;  however,  let  us  show  more  clearly  why  it  is 
an  acceptable  approximation  of  a  true  orbit  in  the  system. 

[  Ref  .  9  :  p  .  2  5 1  ]  The  sequence  of  computed  values,  w  ,  from,  the 
Hencn  map  for  some  initial  condition  y  is  actually  the 
floating  point  representation  of  f(w  ._)  or. 


w-  =  y- 
W  =  f  i  (  f  (V;  )  ) 

w_  =  f  1  (  f  (w  )  ) 
w,  =  f 1 ( f  (w  ) ) 

w, .  ■  =  f  L  (  f  ( w.  )  )  . 

Recall  from  our  definition  that  f(y;)  is  within  e  of  y  .  .  In 
order  for  us  to  show  that  a  computed  orbit,  w. ,  can  be 
arbitrarily  close  to  an  exact  orbit,  x.,  which  5-shadows  the 
orbit,  y  ,  we  must  first  make  an  important  assumption.  We 
assume  that  at  each  iteration  the  error  between  the  computed 
value  of  fly.)  and  the  exact  value  of  fly.)  is  bounded.  That 
is , 

dlw.,fiy  .;))  <  e-  such  that  £■  <  e. 

where  is  some  function  of  the  machine  precision  which 
bounds  the  error  for  all  iterations.  Figure  10  shows  a 
representative  configuration  of  the  orbits  involved.  By  the 
Shadowing  Lemma  dly  ,x  )  <  6. 
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Using  the  triangle  inequality  it  follows  that 

d(w,,x.)  ^  d (Xi,  f  (y,.i)  )  d{Wj,  f  (y... )  )  . 


But , 


d(x, ,  f  (y..,)  )  5  d{x.,y,)  +  d(y. ,  f  (y^.;)  )  . 


Therefore, 

d(w,,x,)<  d(x,,y,)  +  d(y. ,  f  (y,,; )  )  +  d(w_,f;y.;)) 

or 

d(w,,xj^  5  +  e  +  e„. 

If  we  let  5'  =  5  +  e  +  e„  then  it  follows  that,  for  every  6'>e, 
we  can  calculate  a  computed  orbit,  w. ,  that  is  5' -shadowed  by 
an  exact  orbit,  x, .  In  this  way  we  see  that  the  existence  of 
our  computed  orbit  in  the  5 '-shadow  of  an  exact  orbit  actually 
depends  on  the  precision  of  the  machine  in  use. 
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According  to  a  theorem  in  a  recent  paper  (Ref. 15 j, 
specific  criteria  have  been  developed  under  which  we  are 
assured  that  in  the  shadow  of  a  computed  orbit  there  will 
always  exist  an  exact  orbit  within  a  calculable  tolerance. 
The  authors  accepted  that  the  Henon  map  and  the  Henon 
attractor  satisfied  Definitions  1  and  2  previously  mentioned. 
The  accuracy  of  the  numerically  computed  orbits  which  the 
theorem,  measures  applies  to  non-hyperbolic  sets  as  well  as 
hyperbolic  sets  so  compliance  to  Definition  3  may  be  relaxed 
[Ref. 17],  Computations  applying  the  theoremi  in  reference  15 
were  performed  in  Microsoft  Quickbasic  using  a  double 
precision  IBM  system  computer.  The  value  of  e  (the  distance 
between  y  and  f{y,.-. )  which  defines  the  pseudo-orbit)  was  held 
to  approximately  C’’-.  The  authors  applied  the  theorem  and  it 
was  discovered  that  even  after  372,000  iterations  there  is  a 
true  orbit  which  differs  by  at  most  2'-'  (approximately 
1.8626x10''')  from  the  computed  orbit  generated  by  the  Henon 

miap . 

Based  on  the  dense  structure  and  infinite  detail  of  the 
Henon  attractor  it  is  intuitively  reasonable  that  in  the  fcct- 
print  of  the  computed  orbit  there  exists  an  exact  orbit. 
Thus,  the  Shadowing  Lemma  assures  us  that  the  statistical 
evidence  measured  under  computer  analysis  is  significant . 
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II.  PSBTTDORANSCaoaSS  FROM  CHAOS 

A.  LINKINO  CHAOS  TO  THB  BIMARY  DOMAIN 

Symbolic  dynamics  is  a  technique  that  can  be  used  to 
relate  the  dynamics  of  a  particular  system  on  a  metric  space 
to  another  system  on  symbol  space.  Symbolic  dynamics  allows 
us  to  analyze  a  system  by  studying  its  effect  on  symbol  space. 
The  principal  idea  is  that  if  two  systems  are  topologically 
conjugate,  then  their  dynamics  are  •qulvalent. 

DEFINITION:  Suppose  we  have  two  maps:  f:  U  — »  U  and  g:  V V . 
The  functions  f  and  g  are  topologically  conjugate  if  there 
exists  a  homeomorphism  h:U  -4  V  such  that  h  o  f  =  g  o  h.  The 
function  h  is  a  homeomorphism  if: 

i.  h  is  one-tc-one 

DEFINITION:  Suppose  f  is  a  function  from  U  to  V;  f  is 
one-to-ono  if  and  only  if  for  all  elements  U;  and  u^eU, 
f(u.)  =  f{u-)  implies  U;  =  u, . 

ii.  h  is  onto 

DEFINITION:  Suppose  f  is  a  function  from  U  to  V;  f  is 
onto  if  and  only  if  for  any  element  v  in  V  there 
exists  an  element  u  in  U  such  that  v  =  f (u) . 

iii.  h  is  continuous 

iv.  h'‘  is  continuous. 
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Criterion  iv  is  actually  implied  by  criteria  i,  ii,  and  iii. 
Furthermore,  if  f  is  topologically  conjugate  to  g  then  it 
follows  that  g  is  topologically  conjugate  to  f. 

Finding  a  homeomorphism  that  forms  a  topological  conjugacy 
between  two  maps  can  be  extremely  difficult.  However,  this 
does  not  diminish  the  power  of  the  relation.  A  topological 
con:jugacy,  as  we  mentioned,  relates  the  dynamics  of  two 
systems  completely.  [Ref . 20 :p .  27]  For  example,  if  the 
function  f,  from  our  definition,  has  a  period-two  cycle  {p,q}, 
then  ^h'p',h!q)}  is  a  period-two  c:,’cle  for  the  function  g  Tn 
this  way,  all  orbits  for  the  function  f  have  corresponding 
crbits  for  the  function  g.  Furthermore,  if  f  has  a  dense  set 
of  periodic  points  in  U,  then  the  same  is  true  for  g  in  V. 
[Ref . 5 : p . 571 ]  Following  this  reasoning,  if  it  can  be  shown 
that  the  dynamics  of  map  f  exhibit  the  three  characteristics 
of  chaos  previously  mentioned,  and  there  exists  a  topological 
conjugacy  between  f  and  g,  then  map  g  will  also  be  chaotic. 
[Ref .20 :p. 28] 

Can  we  harness  this  chaotic  energy  that  exists  in  one 
system  and  through  an  effective  symbolic  dynamics  transfer  it 
to  another  under  a  homeomorphism  to  achieve  a  useful  result? 
This  is  precisely  the  issue  raised  in  Rejane  Forre's  treatise 
of  November  1990.  [Ref. 12]  Forre  hoped  to  apply  a  presumably 
chaotic  discrete  dynamical  system  to  the  field  of 
cryptography.  Specifically,  she  attempted  to  devise  a  scheme 
to  generate  nearly  random  sequences  of  zeros  and  ones  which 
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could  be  used  for  coding  purposes.  (Any  numbers  generated  by 
a  deterministic  rule  cannot  be  truly  random  and  are, 
therefore,  termed  pseudorandom.)  Forre's  proposed  symbolic 
dynamics  relates  the  apparent  chaotic  dynamics  of  an  orbit  on 
the  Henon  map  to  binary  codespace  by 

h ;  H.  — +  Xc  . 

Here,  A  represents  the  attractor  associated  with  the  mapping 
f  (see  definition  Chapter  I.  section  A.  p.l);  in  partici’la’', 
f  is  the  Henon  map.  The  set  Ss.  also  called  symbol  space  or 
codespace,  represents  the  collection  of  all  infinite  sequences 
of  zeros  and  ones.  Any  chaotic  behavior  exhibited  in  the 
Henon  map  transferred  to  X2  would  be  observed  as  a 
pseudorandom  stream  of  zeros  and  ones  in  L,.  Forre's  symbolic 
dynamics  is  based  solely  on  the  horizontal  component,  x. ,  of 
the  iterates  of  the  Henon  map.  The  elements  of  such  a  binary 
sequence  {S, in  L-  are  defined  as  follows: 

if  Xj  <  x„r-,  then  S.  =  0; 

if  X,  >  x^gr.,  then  S-  =  1, 

for  a  string  of  length  N  where  x„e-  represents  a  dynamic 
median.  (The  dynamic  median  ensures  that  on  the  average  the 
trajectory  will  fall  on  each  side  of  the  median  for  half  the 
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iterates;  this  median  has  been  more  carefully  calculated  as 
Xme-=.4098  in  reference  13).  In  order  to  clarify  the 
computational  process  involved  in  Forre's  symbolic  dynamics 
the  following  example  is  furnished.  Let  us  arbitrarily  choose 
a  particular  (x^.y.)  value  from  the  trapping  region,  say 
(X  ,y  )  =  (1,0)  .  Since  x-  >  x„£.  the  binary  sequence  element 
corresponding  to  is  S,  =  1.  Using  equations  (3.1)  and 
(3.2)  below  we  calculate  ix-.y)  as  follows: 


x„.;  =  1  -  1.4  X,-  +  y, 
y..i  =  •3x„ 


X,  =l-1.4(i)-+0=-.4 
y.  =  .3(1)  =  .3 


Because  x,  <  x^jn  the  next  binary  sequence  element, 
corresponding  to  (x,  ,y,)  ,  is  S;  =  0.  In  this  way  we  can 
calculate  the  entire  forward  sequence  of  binary  elements: 


S  =  {So,S;,S2,S,.  .  .}  =  {1, 0, 1, 0,  .  .  .  }  . 


Although  Forre's  calculated  median  was  inaccurate  (she 
used  a  dynamic  median  x„e..  =  .39912)  her  results  showed  that 
under  the  cryptographic  properties  of  linear  complexity  and 
jump  complexity  the  bitstreams  were  wholly  indiscernible  from 
truly  random  sequences.  However,  she  concluded  that  a  third 
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property,  the  n-tuple  distribution,  was  inconsistent  with 
truly  random  sequences  and  that  this  rendered  the  unaltered 
scheme  unsuitable  as  a  pseudorandom  number  generator. 

B.  RBCBNT  STUDY  OP  FORRB'S  SYMBOLIC  DYNAMICS 

Forre  drew  attention  to  the  poor  n-tuple  characteristics 
of  her  unaltered  symbolic  dynamics.  In  a  more  recent  study 
[Ref.  13]  Heyman  refined  the  calculation  of  the  median,  , 
and  rigorously  investigated  the  n-tuple  property,  as  well  as 
three  other  cryptographic  properties,  to  evaluate  the 
pseudorandom  number  generator.  It  was  concluded  that  the  n- 
tuple  or  runs  property  was  a  minor  detractor  and  that  the 
proposed  use  of  the  scheme  as  a  cryptographic  pseudorandom 
number  generator  was  sound  and  effective .  The  runs  property 
was  mentioned  by  both  Forre  and  Heyman,  but  the  authors' 
conclusions  with  regard  to  the  significance  of  the  property 
were  markedly  different.  Under  these  mixed  results,  the  runs 
anomaly  demands  further  investigation. 

Due  to  the  disparate  conclusions  drawn  by  Forre  and 
Heyman,  it  is  a  natural  next  step  to  evaluate  the  symbolic 
dynamics  under  more  mathematically  rigorous  criteria  in  order 
to  determine  more  objectively  and  conclusively  whether  the 
scheme  is  an  effective  pseudorandom  number  generator.  Because 
the  chaotic  behavior  on  the  attractor  and  certainly  the 
symbolic  dynamics  scheme  itself  are  mathematical  concepts,  we 
evaluate  the  pseudorandom  number  generator  mathematically. 


This  is  certainly  necessary  since  the  previous  papers  were  not 
analytical  and  gave  only  a  cursory  mention  of  the  theory  of 
chaos.  Primarily,  our  task  is  to  collect  evidence  to  prove  or 
disprove  that  this  scheme, 

h  : 

gives  a  homeomorphism.  Our  secondary  objective  is  to  explain 
completely  the  runs  anomaly  in  mathematical  terms.  Thus, 
based  on  the  definition  of  a  homeomorphism  we  will  provide 
experimental  evidence  to  support  the  one-to-one  property. 
Furthermore,  we  will  provide  analytic  proof  that  the  presumed 
homeomorphism  is  not  onto  for  the  proposed  pseudorandom  number 
generator . 


III.  BXPIRIIISIITAL  IVIDINCX  POR  ONX-TO-ONB 

A.  A  CCaCPUTIR  PARASKai 

The  Henon  attractor  A  =  is  an  infinite  set  of  points 
produced  by  the  Henon  map.  The  initial  condition  of  the  Henon 
recurrence  dictates  which  orbit  on  H^  we  produce.  Let  X  ,'K  be 
two  infinitely  long  sequences  of  points  (orbits)  based  on  two 
different  initial  conditions  (s-.t.)  and  (l.,ra-)  respectively, 
where  (s-,t  ),  (l,,rrt)  are  from  H^.  The  sequence  X  will  be 
defined  as  {(s-,t-),  {Sj,t.),  (Sj,  tj)  ,  .  .  . }  .  The  sequence  will 
be  {  f  1- ,  ra.)  ,( 1; ,  md  ,( 1;, m, .  Although  the  two-dimensional 
plots  of  X;  and  X-  are  indistinguishable,  the  two  sequences  of 
points  X;  X;  are  not  the  same.  Let  h(X._),h(Xd  belong  to  Z, 
the  space  of  all  infinite  sequences  of  zeros  and  ones.  For 
the  symbolic  dynamics  h:HA  ^  to  be  one-to-one,  it  must  be 
true  that  if  h(Xd  =  h(X;)  then  X-  =  X;.  That  is,  if  two 
binary  sequences  belonging  to  are  identical  then  they  must 
be  mapped  from  the  same  initial  point  (x,y)  in  H. . 

We  cannot  hope  to  find  with  any  precision  the  points 
belonging  to  H^.  Furthermore,  by  increasing  the  number  of 
possible  initial  conditions  or  keys  from  which  we  can  generate 
orbits  close  to  the  attractor  we  enhance  the  cryptographic 
qualities  of  the  pseudorandom  number  generator.  That  is,  v;e 
reduce  substantially  the  possibility  of  the  key  being  found 
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and  the  code  being  broken.  For  these  reasons  we  modify  our 
symbolic  dynamics  scheme  to  h:TR— That  is,  we  increase 
the  number  of  possible  keys  to  all  those  possible  points  from 
the  trapping  region.  Under  the  Henon  map  the  trapping  region 
is  also  invariant;  therefore,  h:TR— is  a  symbolic  dynamics 
analogous  to 

If  two  binary  sequences  are  not  the  same  then  they  m.ust 
have  originated  from  different  initial  (x,y)  values.  The 
total  count  of  computer  representable  numbers  in  the  trapping 
region,  although  extremely  large,  is  finite.  Thus  the  task  cf 
proving  that  the  homeomorphism  is  one-to-one  on  an  uncountably 
infinite  number  of  points  in  the  plane  is  avoided  and 
unnecessary.  Certainly,  no  irrational  numbers  are  computer 
representable.  Since  all  pseudorandom  number  generators  are 
implemented  in  a  computer  environment  we  qualify  our  goal  to 
showing  that  the  presumed  homeomorphism  is  computer  cne-to-cr.e 
on  TR.  That  is,  we  attempt  to  show  that  the  homeomorphism  is 
one-to-one  with  respect  to  the  computer  representable  points 
in  the  quadrilateral. 

We  propose  then  a  computer  paradigm  that  models  the  finite 
numlser  of  computer  representable  numbers.  The  model  will  also 
give  us  insight  into  the  complicated  dynamics  of  the  map. 
Consider  the  trapping  region  or  quadrilateral  as  a  grid 
consisting  of  a  large  but  finite  number  of  points.  Figure  11 
shows  the  quadrilateral  and  the  median  Xy^-  =  .4098. 
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Flgur*  11  Trapping  region  and  dynamic  median 


•  • 


The  quadrilateral  is  divided  into  two  sides  by  the  median. 
An  initial  point  chosen  from  the  left  quadrilateral  will 
correspond  to  a  sequence  under  hiTR-^Sj  that  will  begin  with 
the  binary  digit  zero  (0).  Similarly,  an  initial  point  from 
the  right  quadrilateral  will  correspond  to  a  sequence  that 
will  begin  with  the  binary  digit  one  (1) .  Let  us  focus  on  the 
right  quadrilateral.  {Our  argument  will  extend  to  the  left 
quadrilateral . ) 

We  model  the  computer  points  in  the  right  quadrilateral 
region  by  first  enclosing  it  in  a  rectangle  of  minimum  area. 
By  choosing  a  particular  spacing  we  can  fill  the  rectangle 
with  equally  spaced  points  to  produce  a  Cartesian  coordinate 
system.  If  we  disregard  the  points  in  the  rectangle  but 
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outs-de  the  trapping  region  we  have  a  grid  which  can  be  used 
to  model  the  large  but  finite  number  of  computer  representable 
points.  If  we  can  show  conclusively  that  as  our  spacing 
becomes  smaller  and  smaller,  the  binary  sequences  of  these 
points  are  all  different  from  that  of  a  point  chosen  at  random 
not  belonging  to  the  grid  points,  then  we  should  be  convincea 
that  h  is  one-to-one.  Fundamentally,  we  first  choose  a  fine 
grid  spacing  for  the  right  quadrilateral  that  corresponds  to 
a  large  number  of  points  in  the  plane.  As  we  iterate  all  the 
grid  points  and  our  random  point  through  the  Henon  map  we  cull 
out  those  points  that  at  each  iterate  do  not  give  the  same 
binary  element  from  {0,1}  as  the  random  point.  We  must 
convince  ourselves  that  there  exists  some  iterate  for  which 
ail  the  binary  sequences  corresponding  to  the  grid  points 
differ  from  the  binary  sequence  of  the  random  point.  Thus,  we 
focus  on  a  succession  of  subsets  of  our  original  grid  pjints 
at  each  iteration.  In  this  systematic  way  with  computer 
assistance  we  hope  to  reveal  the  answer  to  our  question. 

B.  PUTTING  THB  MODEL  TO  WORK 

Let  us  begin  with  an  example  of  our  model  in  practice.  This 
example  was  chosen  because  it  allows  the  reader  a  clear, 
typical  depiction  of  the  procedure  in  a  small  number  of 
Iterations.  Various  MATLAB  programs  for  this  process  on  both 
sides  of  the  quadrilateral  can  be  found  in  Appendix  A  under 
the  names  GRDCOMPl.M  through  GRDCOMP4.M.  As  shown  in  Figure 
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12  che  model  lets  us  fill  the  right  side  of  the  quadrilateral 
with  ecjually  spaced  points  determined  ty  a  specified  value  for 
the  spacing. 


Figure  12  Typical  run  right  quad  model:  before  iteration 


Here  we  use  a  relatively  coarse  spacing  of  0.06.  The  point 
(W,Z)  represented  by  a  cross  within  the  field  is  chosen  at 
random  but  should  be  different  from  the  grid  points.  Figures 
13  through  21  represent  subsequent  iterates  of  the  Henon  map 
where  only  those  points  from  the  grid  which  give  the  same 
binary  sequence  as  the  randomly  chosen  point  (W,Z)  are 
preserved.  The  median,  =  .4098,  is  the  dynamic  median,  as 
previously  mentioned;  therefore,  it  is  not  surprising  that  we 
lose  approximately  half  the  points  at  each  iteration. 


•  • 
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MaximM  of  Poiaa  in  Cnd  wiife  S«aa  Bimiy  Soquoacc 


Figure  13  Typical  run  ric;ht  quad  model:  Iterate  1 


Ma«e<nemofPoiiusuiCinl«iitiSainn  Binoty  S«jooice 


Figure  14  Typical  run  right  quad  model :  Iterate  2 
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Miwiimm  of  Mm  i*  Ghd  with  S«m  liMy  Siyrr 


Figure  15  Typical  run  right  quad  model:  Iterate  3 


Mo»tm«ai  of  Mna  ia  Ghd  wiib  Some  Biaaiy  Sogumer 


Figure  16  Typical  run  right  quad  model:  Iterate  4 
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Figura  17  Typical  run  right  quad  model:  Iterate  5 


Figure  18  Typical  run  right  quad  model:  Iterate  6 
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rlgur«  20  Typical  run  right  quad  model:  Iterate  8 


of  PoiMiia  Cfitf  S«iK  Biaaty  Saqancc 


Figure  21  Typical  run  right  quad  model:  Iterate  9 


By  iterate  10  (not  shown)  no  grid  point  is  on  the  same  side  of 
Xheo  as  (W,Z)  .  Figure  22  gives  us  an  idea  of  how  quickly  the 
points  that  follow  the  orbit  of  (W,Z)  diminish  in  nxrmber.  The 
curve  (1/2)'’  is  plotted  as  a  coir?)arison .  TABLE  1  shows  the 
count  of  points  for  iterates  1  through  9  that  give  the  same 
binary  sec[uence  as  (W,Z). 


TABLE  1  NUMBER  OF  POINTS  THAT  GIVE  THE  SAME  BINARY  SEQUENCE 


NUMBER  OF  POINTS  THAT  GIVE  THE  SAME  BINARY  SEQUENCE 

ITERATE 

0 

B 

2 

3 

B 

5 

6 

B 

8 

9 

NO.  PTS 

98 

65 

13 

9 

B 

B 

4 

B 

B 

B 
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Figure  22  Grid  points  with  same  sequence  after  n 
Iterations 


Presumably  based  on  the  map's  sensitive  dependence  on  initial  •  0 

conditions  there  occurs  a  natural  spreading  effect,  most 

noticeable  in  the  left  half  of  the  quadrilateral,  that  forces 

the  previously  close  points  away  from  each  other.  A  typical  • 

example  of  a  run  of  this  procedure  using  the  left 

quadrilateral  can  be  found  in  Appendix  B. 

Does  there  always  exist  an  iterate  N  where  the  sequence  • 

for  {W,Z)  diverges  from  all  sequences  corresponding  to  even  an 
extremely  fine  grid?  Certainly  this  example  only  gives  a  . 
taste  of  the  capability  of  the  computer  to  fill  the  right  or  • 

left  quadrilateral  with  points.  Based  on  the  limitation  of 
computer  memory  it  is  not  possible  to  use  a  Cartesian  spacing 
small  enough  to  model  even  a  small  fraction  of  the  huge  number  • 
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of  computer  representable  points.  This  is  not  our  goal.  We 
concern  ourselves,  instead,  with  the  trend  in  the  iterate  N 
where  all  the  grid  points  diverge  from  the  iterate  of  (W,Z). 

Let  N  be  the  iterate  at  which  all  the  sequences 
corresponding  to  the  grid  points  diverge  from  the  sequence 
corresponding  to  (W,Z).  Using  our  last  example  we  plot  the 
grid  point  that  lasted  until  iterate  N  (call  it  (u,v))  with 
(W,Z)  for  the  9  iterations.  The  trajectory  of  the  two  points 
is  depicted  in  Figures  23  through  32.  Point  (W,Z),  again,  is 
depicted  as  the  cross  {+)  and  the  point  (u,v)  by  an  (x) . 


Figure  23  (W,Z)  and  (u,v)  prior  to  iteration 


Mavwni  of  K  Qmw  M  HA  aad  Macfeiaf  SaqaMct  ?t  jmI  Oiv«rft 


Figure  24  (W,Z)  and  (u,v) :  Iterate  1 


i 


Figure  26  (W,2)  and  (u,v):  Iterate  3 


Movcmau  of  Pi  Cbowa  aa  HA  Md  Masiuai  Saqnaec  Pt  onifl  OivfT|c 


Figure  27  (W,Z)  and  (u,v) :  Iterate  4  (Notice  that  the 
two  points  are  essentially  superimposed) . 


Figure  30  (W,Z)  and  (u,v):  Iterate  7  (Notice  that  the 
two  points  are  essentially  superimposed) . 


Movodcoi  of  Pt  CheMHi  M  HA  lad  MackiNI  Pt  MMil  Diverte 


Figure  31  (W,2)  and  (u,v) :  Iterate  8 


Notice  Che  proximity  of  (W,Zi  and  (u,v) .  They  are  in  no  ^ 

way  superimposed.  Indeed,  in  a  Euclidean  sense  the  two  points 
are  not  close  at  all.  These  positions  may  seem  acceptable  but 
even  when  the  spacing  between  points  is  dropped  to  the  limits  0 

of  computer  memory  the  proximity  of  (W,Z)  to  (u,v)  does  net 
appreciably  change  from  this  typical  example.  It  is  a  fact 
chat  throughout  the  trajectories  of  the  two  points  there  may  • 

be  iterations  where  the  points  are  virtually  superimposed  as 
in  this  example,  but  as  a  rule  (W,Z)  and  (u,v)  are  not 
initially  neighbors .  Furthermore,  as  in  the  typical  run  using  • 

Che  left  quadrilateral  in  Appendix  B,  sometimes  more  than  a 
single  point  (u,v)  gives  the  same  binary  sequence  as  (W,Z) . 

Typically  one,  but  up  to  seven  points  (u,v),  in  arbitrary/  • 
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observations,  have  mirrored  the  bitstream  of  (W,Z)  and  these 
multiple  matches  likewise  are  not  tightly  packed  points.  This 
phenomenon  lends  support  to  our  model  because  it  is  not  the 
points  closest  to  (W,Z)  which  typically  give  the  same  binary 
sequence.  Typically  many  other  points  were  much  closer  in  a 
Euclidean  sense  to  (W,Z)  but  they  were  rejected.  That  is,  we 
cannot  look  only  in  the  immediate  vicinity  of  (W,Z)  for  points 
that  offer  good  bitstream  matches.  Certainly,  if  a 
cryptographer  knew  that  his  chances  for  regenerating  an 
identical  binary  sequence  were  much  higher  in  the  vicinicy  of 
(W,Z)  he  could  use  it  to  his  advantage.  Therefore,  this 
quality  is  benericial. 

The  GRDCOMPl.M  -  GRCC0MP4.M  models  (mentioned  previously; 
perform  a  grid  comparison  between  the  chosen  point  {W,Z)  and 
the  grid  developed  on  a  particular  side  of  the  dynamic  median. 
The  models  GRDCOMPl.M  and  GRDC0MP3.M  examine  behavior  in  t.he 
right  quadrilateral,  and  GRDC0MP2.M  and  GRDC0MP4.M  model  the 
left  quadrilateral.  The  GRDCOMPl.M  and  GRDC0MP2.M  models  use 
a  predetermined  (W,Z)  value  from  their  respective  side  of  the 
quacrilateral  but  not  on  the  attractor  (H^^)  to  compare  against 
the  developed  grid.  The  GRDC0MP3.M  and  GRDC0MP4.M  models 
differ  from  GRDCOMPl.M  and  GRDC0MP2.M  in  that  (W,Z'  is 
iterated  through  the  Henon  miap  to  ensure  that  the  point 
neighbors  the  attractor  before  the  comparisons  begin. 
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These  computer  models  give  us  a  feel  for  the  N  value  to 


expect  whe:.  we  designate  a  particular  spacing  in  the  right  or 
left  grid  field.  The  computer  models  PROOFl.M  -  PR00F4.M  are 
identical  to  the  GRDCOMPl.M  -  GRDC0MP4.M  models  except  that 
they  allow  us  to  vary  the  spacing.  We  are  able  to  desigr-te 
a  coarse  and  fine  spacing  range  and  an  increment  to  use 
between  these  values.  As  the  Cartesian  spacing  is  steadi'-, 
decreased  and  more  points  fill  the  particular  half -field,  by 
observation,  N  likewise  seems  to  increase  steadily.  This 
behavior  is  pictured  in  Figure  33.  Here,  we  use  the  PROOF^-.M 
model  incrementing  the  spacing  from  .003  to  .5  by  .CCl  for  a 
tota.  of  nearly  500  entries. 
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Although  there  are  a  few  outliera  all  the  N  values  appear  to 
be  clustered  and  to  steadily  increase  overall  as  the  spacing 
becomes  finer  and  finer.  This  behavior  was  observed  for  the 
FR00F2.M  -  FR00F4.M  models  as  well. 

Let  us  establish  .10  as  the  coarse  spacing  upper  bound 
because  it  is  a  round  number  and  it  allows  us  to  fill  either 
half -field  with  more  than  just  a  few  points.  In  computer  runs 
using  this  upper  bound  and  steadily  finer  spacing  ito  the 
limuts  of  computer  memory)  we  observe  the  same  general 
relationship  between  the  decreasing  spacing  and  N  values. 
That  13,  as  the  spacing  is  steadily  decreased  the  value  of  N 
correspondingly  increases. 

To  further  test  our  models  certainly  we  must  not  limit 
ourselves  to  a  comparison  of  a  single  (W,Z)  value  over  the 
range  of  spacing  values.  Because  we  observe  the  same  general 
behavior  in  our  GRDCOMP  models  whether  we  use  a  particular 
iW,Z)  not  on  the  attractor  or  a  (W,Z)  from  the  attractor  let 
us  choose  an  arbitrary  number  of  {W,Z)  values  and  repeat  our 
models  'FRFl.M  -  PRF4.M) .  For  simplicity  the  arbitrar^^  number 
of  (W,Z)  values  (31)  are  taken  from  a  line  belonging  to  the 
particular  side  of  the  quadrilateral.  Figure  34  shows  a  run 
of  model  FRFl.M  (which  corresponds  to  our  GRDCOMPl.M  model  but 
for  multiple  spacings  and  {W,Z)  values). 
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Figure  34  Reduction  in  spacing  vs.  N  (PRFl.M  model) 


This  run  again  is  typical  of  all  runs  of  the  models  PRFl.M  - 
PRF4.M.  This  run  increments  the  spacing  from  .01  to  .10  in 
steps  of  .001.  Notice  that  there  do  not  appear  to  be  31 
values  for  each  spacing  value.  This  phenomenon  actually 
supports  our  observations  because  many  of  the  points  are 
superimposed  on  each  other  on  the  graph.  Again,  we  observe  an 
apparently  steady  increase  in  N  values  as  the  spacing  becomes 
finer  with  few  significant  outliers.  This  graph  contains 
almost  3000  points. 

In  order  to  show  experimentally  that  the  observed 
behavior  exists,  we  select  a  modest  fine  spacing  as  a  lower 
bound  and  a  tiny  increment.  By  using  an  extremely  small 
spacing  we  are  able  to  collect  an  enormous  amount  of  data  to 
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use  as  statistical  evidence.  The  limitation  is  not  con^iuter 
memory  but  the  mainfraune  graphical  statistical  system  (AGSS) 
which  allowed  a  sanple  size  of  45136  points.  For  each  of  the 
31  (W,Z)  values  1456  separate  runs  were  completed.  This 
corresponded  to  an  increment  of  .0000625  between  .01  and  .10. 
The  results  are  pictured  in  Figures  35-38. 
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NONLMeAR  CURVE  FIT 
SCATTER  PLOT.  SS2-45I36 


Flgura  36  Nonlinear  curve  fit  model  PRF2.M 


NONLINEAR  CURVE  FIT 
SCATTER  PLOT.  SS2«45136 


Figure  37  Nonlinear  curve  fit  model  PRF3.M 
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NONLMCAR  CURVE  FIT 
SCAHCR  PLOT.  SS2-4SI36 


Figure  38  Nonlinear  curve  fit  model  PRF4.M 

The  figures  correspond  to  the  models  PRFl.M  to  PRF4.M.  The 
nonlinear  curve  fit  takes  into  account  all  45136  points 
although  only  a  representative  handful  are  pictured.  The 
curves  respectively  correspond  to  the  following  nonlinear 
equations : 


PRFl.M: 

curve 

1 

21.092 

e 

-10.5^3x 

PRF2 .M: 

curve 

2 

24.989 

e 

-3.0599X 

PRF3 .M: 

curve 

3 

18.506 

e 

-U  .692x 

PRF4 .M: 

curve 

4 

20.824 

e 

-9.3852X 

All  four  nonlinear  curve  fits  are  monotonically  increasing 
from  right  to  left.  This  result  seems  to  suggest  that  in  a 
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statistically  acceptable  fashion,  given  a  particular  spacing 
for  our  grid  we  can  predict  an  N  value  for  which  the  binary 
sequence  of  (W,Z)  will  diverge  from  those  of  the  grid  field. 

How  small  can  we  make  c  r  grid  spacing  before  our  models 
fail?  Certainly,  for  any  particular  computer  there  will  exist 
a  value,  call  it  e,  such  that  computer  numbers  within  that 
tolerance  will  be  considered  the  same.  That  is,  if  x  and  y 
are  computer  representable  numbers  and  |  x-y  |  <  e,  then  x  =  y. 

The  value  of  £  varies.  Let  computer  epsilon  be  defined  as 
a  value,  e,  below  which  1  +  e  =  1;  386  Matlab  version  3 . 5M 
computes  this  value  as  approximately  2.24  x  10' '.  However, 
the  value  of  e  (using  the  same  definition)  such  that  0  +  e  = 
0  is  on  the  order  of  10'’°®  for  the  same  system.  One  can 
speculate  that  at  some  diminutive  spacing  S,  our  model  will 
fail.  That  is,  there  will  exist  a  spacing  where  the  computer 
will  not  be  able  to  build  the  grid.  The  computer  model  will 
calculate  the  first  point  x  in  the  rectangle  which  encompasses 
the  quadrilateral.  But,  because  the  spacing  S  is  smaller  than 
the  £  tolerance  that  the  computer  needs  to  recognize  the  next 
point  y,  X  and  y  are  seen  as  the  same  point.  Only  the  first 
grid  point  will  appear  in  the  grid.  However,  this  is  in 
keeping  w'ith  cur  definition  of  one-to-one  for  our  presumed 
homeomorphism  because  if  two  sequences  are  wholly  the  same 
then  they  come  from  the  same  initial  condition.  That  is,  the 
sequences  are  generated  from  two  initial  conditions  (points! 
that  the  computer  recognizes  to  be  the  same  in  finite 
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precision  although  the  initial  conditions  would  differ  in 
infinite  precision.  Thus,  we  have  experimentally  shown  that 
for  the  symbolic  dynamics  if  we  have  two  unequal  initial 
conditions  then  after  some  iteration  N  the  corresponding 
sequences  should  diverge.  Therefore,  we  have  experimental 
evidence  to  support  the  statement  that  the  symbolic  dynamics 
IS  one-to-one. 
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IV.  ONTO 


X.  ONTO  PROPERTY 

Whether  or  not  the  proposed  symbolic  dynamics  has  the  onto 
property  previously  mentioned  is  directly  brought  into 
question  by  the  runs  anomaly  which  was  closely  examined  in 
reference  13.  This  property,  defined  in  Chapter  II.  section 
A.,  informally  states  that  under  the  symbolic  dynamics,  every 
possible  sequence  of  zeros  and  ones  is  possible  or  realizable . 
If  a  pseudorandom  number  generator  has  a  favorable  runs 
property,  all  possible  n-tuples  (sequences  of  zeros  and  ones 
of  length  n)  must  not  only  be  possible  but  their  occurrences 
must  be  balanced.  The  unbalanced  count  of  particular  n-tuples 
is  precisely  the  runs  anomaly  pointed  out  by  Forre  as  the 
fatal  flaw  of  her  pseudorandom  number  generator. 

The  runs  property  of  a  binary  sequence  is  tested  by 
counting  the  occurrences  of  the  2'  different  possible  n- 
tuples.  For  example,  in  Heyman's  study,  typical  binary 
sequences  of  length  10'‘  were  used  to  test  2 -tuples  cn  a 
computer.  There  are  2-=4  2 -tuples:  {0,0}, {0,1}, {1,0}, {1,1}. 
The  program  ONTO.M  (Appendix  A)  counts  the  four  different  2- 
tuples  in  a  binary  sequence.  Figure  39  shows  a  bar  graph 
which  gives  the  count  or  the  four  z -tuples  for  a  typical 
binary  sequence  of  length  5  x  10'. 
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Flgur*  39  Incidence  of  2 -tuples  in  typical  Henon 
generated  binary  sequence 

Notice  that  the  horizontal  axis  is  actually  the  decimal 
equivalent  of  the  n-tuple  shifted  by  one  from  zero.  That  is, 
instead  of  the  following  2 -tuple  to  decimal  correspondence; 


2-tuple 


decimal 


we  use 


00 

01 

10 

11 


1 

2 

3 

4. 


We  will  use  this  shifted  correspondence  for  all  n-tuples.  The 
unbalanced  or  uneven  bins  of  2 -tuples  in  Figure  39 
demonstrates  the  runs  anomaly.  The  depth  of  this  anomaly  is 
shown  in  Figures  40-42  which  correspond  to  3,4  and  5-tuple 


■5; 


00 

- 

0 

where 

[0 

0] 

*  [2 

1]‘  =  0 

• 

• 

01 

- 
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1] 

•  [2 

II 

10 

- 

2 
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[1 

0] 
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- 
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1] 
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counts  respectively  for  typical  binary  sequences  of  length  5 


X  10^.  Modifications  of  ONTO.M  were  used  to  produce  the 
associated  data  for  finding  the  count  of  2-tuples  to  17- 
tuples  . 


Figure  40  Incidence  of  3 -tuples  in  typical  Henon 
generated  binary  sequence 
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Shifted  Declaal  Equivalent 

Figure  42  Incidence  of  5-tuples  in  typical  Henon 
generated  binary  sequence 

Notice  that  in  Figure  41  the  bin  corresponding  to  the  decimal 
13  is  empty.  The  bin  corresponding  to  the  decimal  7  is  also 
empty  but  in  tests  using  binary  sequences  of  length  10“  this 
bin  is  not.  The  13  bin,  however,  remains  empty  in  the  10“ 
length  test. 

It  was  believed  (Ref. 13]  that  as  sequences  of  greater  and 
greater  length  were  tested,  the  missing  sequences  would  be 
found  although  the  runs  property,  namely  the  balance  of  the 
bins,  probably  would  not  improve.  Under  memory  constraints  of 
the  runs  property  test  using  a  Sparc  station  2,  the  length  of 
the  longest  possible  testable  binary  sequence  was  roughly  10“ . 
Using  initial  conditions  from  all  four  quadrants  in  separate 
tests  with  binary  sequences  of  length  10“  no  4-tuple  was  ever 
found  to  correspond  to  the  decimal  13.  The  decimal  13 
corresponds  to  the  binary  4-tuple  {1,1, 0,0).  It  is  therefore 
expected  that,  when  the  5-tuple  runs  are  tested,  those  5- 


tuples  which  contain  the  4-tuple  {1,1, 0,0}  (see  Table  2) 
should  be  empty. 


TABLE  2  5-TDPLBS  CONTAIMIMG  {1,1, 0,0}  AND  DBCIMAL  BQUIVALBNTS 


r 

5-tuples  with  (1,1, 0,0} 

decimal  equivalent 

O 

O 

o 

25 

O 

O 

26 

(C, 1, 1, 0, 0} 

13 

{1,1, 1,0,0} 

29 

This  is  indeed  the  case;  however,  these  are  not  the  only 
subsequences  of  length  5  unrealized  in  a  sequence  of  length 
10'.  Table  3  lists  those  additional  5-tuples  that  are  not 


realized  and  their  corresponding  decimal  values. 


TABLX  3  S'TUPLKS  AMD  THSZR  DKCIMAL  BQaiVALXNT 


5 -tuple 

decimal  equivalent 

{0, 0, 0, 0, 0} 

1 

o 

o 

O 

o 

5 

{1,  0, 1, 1, 0} 

23 

{1- 1, 0, 1, 1} 

28 

In  fact,  the  property  was  tested  for  n-tuples  from  n  =  2 
to  li"  (see  Appendix  C  for  6  through  16-tuple  runs)  and  an 
increasing  percentage  of  unrealized  sequences  occurs.  Table 
4  shows  the  number  of  unrealized  sequences  and  the  total 
number  of  possible  n-tuples  for  a  particular  n. 
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TABLE  4  OHRBALIZBD  SEQUBMCBS  AND  TOTAL  POSSIBLE  n-TUPLES 


n-TUPLE 


UNREALIZED 


SEQUENCES 


POSSIBLE 


PERCENT 


SEQUENCES  UNREALIZEI 


98.36 


Figure  43  shows  how  the  percentage  of  unrealizable  sequences 
increases  as  the  n-tuple  length  increases. 


B.  ANALYTIC  PROOF  OF  UMRSALIZABLB  FOUR-TUPLE 

Based  on  the  previous  data  alone,  the  presumed 
homeomorphism  cannot  be  disqualified  from  having  the  onto 
property.  An  analytically  substantiated  example  of  a  wholly 
unrealizable  sequence,  however,  would  suffice.  Let  us 
investigate  the  first  suspected  unrealizable  sequence, 
{1,1, 0,0}.  From  the  runs  tests  we  see  that  the  sequence 
{1,1,0,!}  is  possible  although  (1,1, 0,0}  is  not.  Figure  44 
shows  the  portion  out  of  5000  points  on  the  attractor  which 
gives  the  sequence  {1,1,0}. 


A 
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Figur«  44  Points  from  a  sequence  of  5000  chat  give  the 
sequence  {1,1,0} 

These  are  also  the  same  points  which  give  the  sequenc 
{1,1,0,!},  since  none  gives  the  sequence  {1,1, 0,0}.  It  i 
apparent  that  these  points  are  localized  in  a  particular  area 
Figures  45-47  show  the  subsequent  iterates  of  these  points. 


Figur*  45  First  iterate  of  points  that  give  the  sequence 

(1,1,0} 


Fir'ire  46  Second  iterate  of  points  that  give  the 
sequence  {1,1,0} 


{1,1,0} 

• 

It  is  clear  that  all  these  points  under  iteration  correspond 

to  the  sequence  {1,1,0,!}.  A  quadrilateral  is  placed  around 

the  points  of  Figure  44  in  Figure  48  (call  it  the  • 

subquadrilateral  or  subquad)  . 


Figur*  48  Subquadrilateral  placed  around  those  points 
that  give  the  sequence  {1,1,0} 


By  iterating  the  points  belonging  to  this  subguad  using  the 
prograun  FINDllO.M  in  Appendix  A  (see  Figures  49-51)  we  see 
computer -genera ted  evidence  that  no  points  within  the  subquad 
correspond  to  the  binary  sequence  (1,1, 0,0}. 


Flgur*  49  First  iterate  of  points  comprising 
subquadrilateral 


Figure  50  Second  iterate  of  points  comprising 
subquadrilateral 


Figure  51  Third  iterate  of  points  comprising 
subquadrilateral 

This  exclusion  zone  does  plainly  suggest  that  the  sequence 
{1,1, 0,0}  is  unrealizable.  Using  a  binary  Henon  sequence  of 
length  lOS  approximately  5115  points  lie  in  the  subquad  or 
exclusion  zone.  This  suggests  that  an  orbit  has  a  probability 
density  of  approximately  5.115%  in  the  subquad.  Since  the 
Henon  map  has  been  shown  numerically  to  have  topological 
transitivity,  we  expect  that  under  reverse  iteration  of  each 
point  in  the  subquad  there  will  be  a  preimage  which  will  lie 
in  the  subquad.  Therefore,  this  means  that  all  points  pass 
through  the  subquad  and  that  no  points  which  pass  through  this 
window  give  { 1 , 1 , 0 , 0 } . 


Let  us  look  more  closely  now  at  the  points  that  would  give 
{1,1, 0,0}  instead  of  those  that  would  not.  If  we  can  show 
that  the  set  of  points  that  can  possibly  give  (1,1, 0,0}  is  the 
empty  set,  then  we  will  have  proved  that  this  proposed 
symbolic  dynamics  does  not  possess  the  onto  property.  Recall 
the  Henon  map  introduced  in  Chapter  I: 

X, .  =  1  -  1 . 4x.  -  +  y. 

y..;  =  .3x. . 

It  follows  that; 

x...  =  1  -  1.4xv  +  .3x..;.  (4.1) 

Our  first  goal  is  to  find  the  solution  set  that  corresponds  to 
{1,1,0}.  In  order  for  a  sequence  of  (x,y)  values  to 
correspond  to  {1,1,0}  the  x  values  {Xr,'.  x. ,  x...}  must  obey  the 
following  inequalities: 


X,.; 

> 

(4.2) 

> 

(4.3) 

X, 

< 

(4.4) 

where  x„^_-=.4098,  the  dynamic  median  rounded  to  four 
significant  digits  (previously  discussed).  Inequality  (4.2; 
can  be  substituted  in  equation  (4.1)  as  follows: 

X...  >  1.12294  -  1.4x/  (4.5) 
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where  the  constant  1.12294  is  rounded  to  five  significant 
digits  as  an  underestimate  in  order  to  make  our  proof  valid. 
Moreover,  in  all  such  cases  we  round  in  the  appropriate  (i.e. 
conservative)  direction.  In  Figure  52  we  see  via  (4.5)  that 
the  solution  set  corresponds  to  the  region  above  the  inverted 
parabola,  in  the  XnXn,; -plane. 
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Figure  53  shows  the  possible  solution  set  under  restrictions 
(4.2), (4.3)  and  (4.4) 


{1,1,0}  under  restrictions 

By  inequality  (4.4),  has  a  maximum  value  less  than  or 
equal  to  x^gj.  If  we  let  x.,..  equal  exactly  x^g^  in  inequality 
(4.5)  then  x„  is  further  restricted  to  a  value  greater  than 
.71371.  Thus,  the  solution  set  corresponding  to  equation 
(4.1)  is  limited  to  the  region  shown  in  Figure  54  and  reviewed 
as  follows: 

X.,.;  :  region  above  the  inverted  parabola  :  1 
X,  :  region  such  that  x.,  >  .71371  :  i 
X...  :  region  such  that  x.,..  £  .4098  :  0 


71 


•  • 


•  t • 


A 


sequence  {1,1,0} 


9) 


♦ 
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By  equation  (4.1)  it  follows  that  if 

x„.;  =  1  -  1.4x„^  +  .3x„.; 

then 

x,.2  =  1  -  1.4x„,i^+  .3x„.  (4.6) 

To  achieve  the  sequence  {1,1, 0,0),  Xn.2  must  be  less  than  or 
equal  to  .4098.  However,  it  is  clear  that  the  largest  x.„,_ 
value  in  equation  (4.6)  occurs  at  a  minimum  value  for  x.,  and 
a  maximum  value  for  Xp.,,  which  lies  in  the  solution  region  in 
Figure  54 .  The  minimum  x„  value  and  maximum  x.,.,  value  are 
respectively  .71371  and  .4098.  Using  these  optimum  values  in 
equation  (4.6) , 
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X...  -  1  +  .3(. 71371)  -  1.4(.4098)- 
-  .979003. 

The  symbol  (-)  indicates  that  we  round  to  the  number  cf 
significant  digits  shown  and  that  we  round  in  the  appropriate 
(conservative!  direction  as  previously  defined.  This  x, . 
value  (not  surprisingly)  corresponds  to  the  point  with  the 
minimum  x  value  in  Figure  51.  that  is,  the  third  iterate  of 
the  subquad  region  which  only  gives  the  sequence  {1,1, 0,1). 
Of  course,  this  minimum  x..^  value  exceeds  .4098,  and 
therefore,  the  sequence  {1,1, 0,0}  is  not  realizable.  Thus,  we 
have  found  a  counterexample  to  the  proposed  homeomorphism. 
Furthermore,  it  is  shown  in  Appendix  D  that  the  sequence 
(0,0, 0,0,0}  IS  also  analytically  unrealizable. 
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V.  CONCLOSIOHS 


A.  RUMS  AM(»fALY-FATAL  FACTOR  OR  LIMITIHO  FACTOR? 

We  have  seen  that  the  classical  Henon  map  exhibits  the 
attributes  normally  accepted  as  characterizing  chaotic 
dynamical  systems.  These  properties  include  sensitive 
dependence  on  initial  conditions,  topological  transitivity  and 
a  dense  set  of  periodic  points.  Despite  this,  however,  we 
have  provided  ample  numerical  evidence  and  rigorous  analytical 
proof  that  the  proposed  symbolic  dynamics  scheme  h:TR— for 
generating  pseudorandom  binary  sequences  is: 

1)  not  a  homeomorphism  for  the  proposed  symbolic  dynamics 
scheme  since  it  is  not  onto.  (Because  the  Henon  attractor  is 
a  subset  of  TR  we  can  also  conclude  chat  h:H;^  — »  I.  is  not 
onto .  ) 

2)  highly  restricted  in  its  viability  as  a  pseudorandom 
n'umber  generator. 

We  have  shown  that  not  all  binary  sequences  are  generated  with 
equal  frequency.  In  particular  we  have  shown  that  certain 
sequences  are  not  realizable  and  that  others  are  very  sparsely 
attained.  These  facts  support  the  observations  and 
conclusions  of  Rejane  Forre  which  suggest  that  this  scheme  is 
unsuitable  as  a  reliable  means  of  generating  pseudorandom, 
numbers . 


T" 


*) 


♦ 


•  • 


74 


•  •  •  •  •  •  •  •  •• 


The  evidence  suggests  that  the  runs  property  severely 
limits  the  potential  of  the  scheme  for  widespread  practical 
use.  It  IS  shown  in  reference  18  that  subsequent  elements  of 
the  binary  sequence  may  actually  be  predicted  with  "'C-8C% 
accuracy  by  an  artificial  neural  network  (compared  to  5C% 
accuracy  for  a  ccin  flip),  most  likely  due  to  the  severity  cf 
the  runs  anomaly.  However,  because  of  its  simplicity  a.nd 
nonlinearity  there  may  exist  som.e  applications  for  which  the 
scheme  would  be  well-suited.  It  is  our  belief  that  the 
general  idea  of  using  a  chaotic  discrete  dynamical  system,  to 
generate  pseudorandom  binary  sequences,  however,  has  merit  and 
deserves  additional  study. 

B.  FUTURE  WORK  FOR  AM  IMPROVED  SYMBOLIC  DYNAMICS 

Our  numerical  results  reflect  the  basic  structure  of  the 
attractor.  Despite  the  lack  of  a  homec  morph  ism.  the  structural 
nature  of  the  attractor  is  apparent  in  the  binary  output  of 
the  symbolic  dynamics  scheme.  Although  the  attractor  is 
accepted  as  possessing  chaotic  attributes,  we  believe  that  due 
to  the  structure  of  this  "chaotic  driver"  (the  classical  Hencn 
map)  pseudo  randomness  is  not  fully  realized  in  this  schem.e. 
We  believe  that  a  sequence  of  more  than  four  zeros  is  not 
possible  because  the  left  fixed  point  is  both  repelling  and 
not  on  the  attractor.  Contrarily,  it  has  been  shewn  that  a 
sequence  of  up  to  23  ones  is  possible  [Ref .16]  .  We  conjecture 
(as  suggested  by  computational  evidence)  that  the  sequences  of 


75 


TT 


9) 


♦ 


•  • 


•  ••••• 


ones  are  possible  because  the  right  quadrilateral  fixed  pcir.': 
13  a  saddle  point.  Specifically,  the  sequence  of  ones  is 
possible  due  to  the  attracting  axis  of  the  saddle  point.  This 
suggests  that  the  symbolic  dynamics  scheme  could  be  effective 
if  there  were  such  a  saddle  point  on  the  attractor  cn  both 
sides  of  the  dynamic  median  . 

Since  the  location  and  nature  of  the  fixed  points  for  the 
Henon  map  depend  on  the  (a,b)  parameters  of  the  m^ap, 
i' Ref  .  1 ; p  .  17 0  ]  there  may  exist  parameter  sets  which  provide 
this  structure.  The  chaotic  bands  present  in  the  bifurcation 
diagrams  suggest  there  exist  many  other  (a,b)  pairs  that  may 
give  even  more  complicated  dynamics.  Proving  the  existence  cf 
a,b'  parameters  that  give  us  these  desired  characteristics  is 
only  the  first  step.  It  is  also  required  that  a  bounded 
attractor  exist  for  the  associated  parameters.  Furthermore, 
Ferre  had  the  luxury  of  being  provided  a  num.erically 
taltulated  trapping  region  from  which  to  take  the  pseudorandom, 
numioer  generator  key.  This  region,  if  it  exists  at  all,  v.'tu.d 
have  to  be  recalculated  for  the  new  ia,b.  pair. 

If  It  could  be  proved  that  there  does  not  e.xist  an  a,b 
pair  which  gives  rise  to  saddle  points  that  possess  the 
previously  m.entioned  attributes,  then  it  m.ay  be  possible  t: 
use  another  map  which  corresponds  to  a  different  attract:! 
with  these  qualities.  However,  the  mater  asset  of  the  .He.non 
recurrence  is  its  utter  simplicity  whic.h  translates  to  the 
fas':  generation  cf  pseudorandomi  sequences. 
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Basing  the  symbolic  dynamics  scheme  on  the  x  value  of  each 
Iterate  is  simple  but  it  is  not  clear  that  it  is  the  most 
effective.  If  we  are  willing  to  accept  a  more  complicated 
scheme  we  could  base  our  split  on  a  different  linear  median  cr 
a  nonlinear  median.  The  goal  in  choosing  a  given  median  would 
be  to  bring  parity  to  the  runs  property  for  the  system  of  dual 
saddle  points  while  retaining  the  other  properties  previously 
mentioned.  A  more  complicated  scheme  might  generate  binar\' 
sequences  more  slowly,  but  the  loss  of  speed  may  be  warranted 
to  improve  the  runs  property. 

It  may  be  argued  chat  despite  an  improved  system  there 
will  still  exist  binary  sequences  that  are  unrealizable. 
However,  the  severity  of  the  associated  runs  anomialy  could  be 
diminished  should  the  shortest  of  these  unrealizable  sequences 
be  of  sufficient  length.  That  is,  we  have  seen  chat  an 
unrealizable  4-tuple  was  catastrophic  in  this  case  because  the 
problem,  of  unrealizability  translated  to  every  subsequent  n- 
tuple.  In  fact,  we  must  anticipate  that  as  described  in 
Chapter  IV  even  more  n-tuples  will  be  unrealizable  than 
expected.  The  degradation  of  the  system  thus  is  directly 
related  to  the  length  (L)  of  the  first  unrealizable 
sequence (s)  .  Since  a  pseudorandom  number  generator  is  yudged 
by  more  than  just  the  runs  property  criterion  (see  conclusion 
Ref. 13)  a  system  which  possesses  a  "large"  L-value  may  still 
provide  an  acceptable  pseudorandom  number  generator. 


APPENDIX  A:  MAIN  PRCX3RAMS 


m 

NOTE:  GRDCOMP  PROGRAMS  REQUIRE  THE  HENREAL  PROGRAM  TO  RUN  Ar 

%  function  [xmat,  ymat]  *  grdconpl (sp, W, Z) 

%  This  function,  takes  a  point  (W,  Z)  from  the  RT  quad,  init  lates 
a  quad  grid  • 

%  field  based  on  a  certain  spacing,  then  iterates  those  grid 
field  points 

%  that  match  the  binary  string  of  (W, Z) [use ( . 1, . 1 )  for  now. 

■ns  long  as 

%  the  strings  match,  points  are  iterated  using  the  Henon 

recurrence.  Only  that  • 

%  (those)  points  that  completely  match  are  finally  plotted 
according  to  the 
%  following  scheme: 

%  b.  original  grid  field  points 

%  r*  original  ( W, Z ) 

%  go  iterated  field  points  with  matching  binary  string  • 

at  nth  Iteration 

%  w+  iterated  (W,Z)  at  nth  iteration 

%  gx  original  xvec,yvec  point (s)  which  matches  binary 

%  string  of  W,Z.  We  also  show  how  W, Z  and 

xvec ( indices ), yvec ( indices )  "walk"  to  the  Henon  attractor  at 

each  Iteration  (with  the  attractor  on  screen) .  •  • 

0=1 inspace (-1.33,1.32,500) ; 
s=-. 1083*0  +  .276; 
u=linspace (1.32,1.245,500) ; 
v=3.64*u  -4.6718; 
g=linspace (-1.06,1.245,500); 

h=.1533*g  -.3344;  • 

e=linspace (-1.06, -1.33,500) ; 
f=-3.407*e  -4.1119; 
split  =  .4098; 

a  =  [-.1083  -l;-3.64  1;-.1562  1]; 
b  =  [-. 2760;-4. 6718;-. 3344] ; 

c  =  [split:3p;1.32] ;  • 

d  =  [-.6:sp: .5] ; 

[x,  y]  =  henreal ' 750 , -1 . 0 , -  .  25 ) ; 

xline  =  [split  split];  yline  =  [-1.32  1.32]; 

lenc  =  length(c);  lend  =  length(d); 

xvec  =  zeros { 1 , lenc*lend) ;  yvec  =  zeros ( 1 , lenc* lend i ; 

if  min ( (a* [w; z] )  >  b)  ==  0  I  w  <  split  • 

disp (' initial  value  is  not  in  the  right  quadr') 
return; 

end 
k  =  0; 

for  p  =  l:lenc 
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for  r  =  l:lend 
k=k-^l; 

if  max ( (a* [c (p) ;d (r ) ] )  <  b)  ==  0 
xvec{k)  =  c(p);  yvec(k)  =  d(r) ; 
else 

xvec(k)  =  -10;  yvec(k)  =  -10; 

end 

end 

end 

xvec  =  reshape (xvec , lenc* lend, 1 ) ; 
yvec  =  reshape (yvec, lenc*lend, 1 ) ; 

m  =  find(xvec  ~=  -10  I  yvec  ~=  -10); 
xvec  =  xvec(m);  yvec  =  yvec(m); 

XX  =  xvec ;  '/y  -  yvec ; 

ww  =  w ; 

zz  =  z; 

axis  1  [-1.32  1.32  -.6  .5]  )  ; 

indices  =  1 : length (xx) ;  %  [1234  .  .  .] 

1  =  1; 

qvec  =  length{xx); 
q  =  indices; 

lenq  =  length(xx)  % 

initialization  of  q 

plot (xvec , yvec , ' b . ' , W, Z, ' r* ' , xline, yline, ' w- ' ) ; hold  on; 
plot ( o, s , ' w. ' , u, V, ' w . ' , e, f , ' w. ' , g, h, ' w. ' ) ;  hold  on; 
title('Grid  field  and  W, Z' ) ;pause;cig;hold  off; 
while  length (q)  >  1; 

axis! [-1.32  1.32  -.6  .5] )  ; 

xO  =  yy  +  ones  ( size  (yy)  )  -  1.4*xx.''2; 

yO  =  .3*xx; 

wO  =  zz  +  1  -  1.4*wv;.''2; 
zO  =  .3*ww; 

if  wO  >  split 

qtemp  =  find(x0  >  split  I  xO  ==  split); 
else 

qtemp  =  find(x0  <  split); 
end; 

if  length (qtemp )  ==  0, break; 

else  q  =  qtemp; 

end 

i  =  i  +  1;  %  counts  iterates  where  at 

least  one  point  matches 
lenq  =  length(q) 

qvec  =  [qvec, length (q) ]; %  plot  of  how  #  with  same  binary- 
sequence 

indices  =  indices (q);  %  decreases  each  time  through  loop 

XX  =  x0(q);  %  preserves  to  next  iterate  those 

matching  values  of 

yy  =  y0(q) ;  %  xO  and  yO 

'WW'  =  w  0  ; 
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zz  =  z 0  ; 
pause; 

plot (XX, yy, 'g. ' ,wO,zO, 'w+' , xliue.yline, ' w- ' ) ; hold  on ; 
plot(o,s, 'w.',u,v, 'w.',g,h, 'w.') ;hold  on; 
title ( 'Movement  of  points  in  grid  with  same  binary 
sequence ' ) 

pause;  clg;hold  off; 

end 
axis  ; 

vec  =  (l/2).''((l:length(qvec))-l); 

plot  (  ( 1  :  length  ( qvec )  )  - 1 ,  qvec,  '  r '  ,  ( 1 ;  length  (qvec )  )  -1 ,  q\'ec  ( 1  >  * 
vec ,'b'); 

t it le f ' Decrease  in  #  of  grid  points  with  same  binary 
sequence  vs  l/2''n') 
pause ; clg ; 

%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  original 
xvec  and  yvec . 

axis( [-1.32  1.32  -.6  .5]); 

plot(xline,yline, 'w',x,y, 'b. ',o,s, 'w. ',u,v, 'w. ',e,f, 'w. ' ,g,h 
, 'w. ' ) ;hold  on; 

plot (xvec ( indices ) , yvec ( indices ) , ' gx' , W, Z, ' r* ' ) ; hold  on; 
title ('W,Z  and  the  grid  point  with  same  binary 
sequence' ) ; pause; 
hold  o f f ; 

xnew  =  xvec ( indices ) ; 

ynew  =  yvec ( indices ) ;  %  now  we  show  a  plot  of  how  W, Z  and 
the  point 

%  that  generates  the  same  binary  sequence 

xnew, ynew 

%  walk  around  the  attractor 
axis ( [-1.32  1.32  -.6  .5]); 

plot(xline,yline, 'w',x,y, 'b.',o,s, 'w.',u,v, 'w. ',e,f, 'w. ' ,g,h 
,  '  w .  '  )  ; 
hold  on; 


xnewt  = 

xnew; 

ynewt  = 

ynew  ; 

wt  = 

w; 

zt  = 

z; 

for  n  = 

l:i  -  1 

xOt 

=  1.0  -  1 . 4*xnewt . ''2  +  ynewt; 

yCt 

=  .3*xnewt; 

wOt 

=  1.0  -  1.4*wt.''2  +  zt; 

zOt 

=  . 3  *wt ; 

wt 

=  wOt; 

zt 

=  z  0 1  ; 

xnewt 

=  xOt; 

ynewt 

=  yOt; 

plot (xOt , yOt , ' gx' , wOt , zOt , ' r* ' ) ; 

♦ 


•  • 
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t It le ( 'Movement  of  pt  chosen  on  HA  and  matching  sequence  pt 
until  diverge' ) ; 

pause;  plot(x0t,y0t, 'ix' ,w0t,z0t, 'i*' ) ; 


%  function  [xnat,  yaat]  ■  grdcoapS (sp/W, Z) 

%  This  function  takes  a  ;W,Z)  from  the  left  quad, init lates  a 
quad  grid  field 

*  based  on  a  certain  spacing,  then  iterates  those  grid  field 
points  that 

%  match  the  binary  string  of  >W,Z)  [use  (-l,-.25)  for  ncwj  . 
As  long  as  the 

%  strings  match,  points  are  iterated  using  the  Henon 
recurrence.  Only  that 

%  (those)  points  that  completely  match  are  finally  plotted 
according  to  the 
%  following  scheme: 

%  b.  original  grid  field  points 

%  r *  original  (W, Z ) 

%  go  Iterated  field  points  with  matching  binary  string 

at  nth  iteration 

%  w+  Iterated  (W, Z)  at  nth  iteration 

%  gx  original  xvec,yvec  point (s)  which  matches  binary' 

string  of  W, Z 

%  We  then  plot  the  iterates  of  W, Z  and  the  point  that  most 

closely  matches 

%  its  binary  sequence. 

o=linspace (-1.33,1.32,500) ; 

s=-.1083*Of .276; 

u=linspace (1.32,1.245,500); 

v=3 .64*u-4.6718; 

g=linspace (-1.06,1.245,500)  ; 

h=.1533’'g-.3344  ; 

e=linspace (-1.06,-1.33,500); 

f=-3 .407»e-4 .1119; 

split  =  .4098; 

a  =  [3.4074  1;  -.1083  -1;  -.1562  1]; 
b  =  [-4.1119;  -.2760;  -.3344]; 

C  =  L -1 . 32 : sp : split ]  ; 

d  =  [-.6:sp:  .5]  ; 

[x,  y]  =  henreal (750 , -1 . 0 , - . 25) ; 

le..c  =  length(c);  lend  =  length(d); 

xline  =  [split  split];  yline  =  [-1.32  1.32]; 

xvec  =  zeros ( 1 , lenc*lend) ;  yvec  =  zeros ( 1 , lenc* lend) ; 
if  min  (  (a*  [’w;  z]  )  >  b)  ==  0  I  w  >  split 

disp ( ' init ial  value  is  not  in  left  quadr') 
return ; 

end 

k  =  G  ; 

for  p  =l;lenc 

for  r  =  l:lend 
k=k+l ; 

if  max ( (a* [c (p ) ; d ( r ) ] )  <  b)  ==  0 
xvec(k)  =  c(p);  yvec(k)  =  d(r); 
else 

xvec{k)  =  -10;  yvec(k)  =  -10; 


TT 


♦ 


•  • 
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t  •  • 


1 _ 


1 


end 


end 

end 

xvec  =  reshape (xvec, lenc*lend, 1 ) ;  %makes  a  column  vector 
yvec  =  reshape (yvec, lenc*lend, 1 ) ; 

m  =  find(xvec~=  -10  I  yvec~=  -10); 
xvec  =  xvec(m);  yvec  =  yvec (m) ; 

XX  =:  xvec ;  yy  =  yvec ; 

WW  =  W/’  ZZ  —  Z/ 

axxs( [-1.32  1.32  -.6  .5] )  ; 

indices  =  1 : length (xx) ;  %  [1  2  3  4  .  .  .] 

1  =  1; 

qvec  =  length(xx); 

q  =  indices;  % 

initialization  of  q 
lenq  =  length (xx) 

plot (xvec , yvec b W, Z r* xline<yline, 'w- ') ;  hold  on; 
plot'o.s, 'w.',u,v, 'w.',e,f, 'w.',g,h, 'w.');hold  off; 
titiel'Grid  field  and  W, Z '); pause; clg;hold  off; 
while  length (q)  >  1; 

axis(  [-1.32  1.32  -.6  .5] )  ; 

xC  =  yy  +  ones  ( size  (yy)  )  -  1.4*xx.''2; 

yO  =  .3*xx; 

wO  =  zz  +  1  -  1.4*ww.''2; 
zO  =  .3*ww; 

if  wO  >  split 

qtemp  =  find(x0  >  split  I  xO  ==  split); 
else 

qtemp  =  find(x0  <  split); 
end; 

if  length ( qtemp)  ==  0, break; 

else  q  =  qtemp; 

end 

1=1+1; 
lenq  =  length(q) 
qvec  =  [qvec, length (q) 1 ; 
indices  =  indices (q); 

XX  =  x0(q);  %  preserves  to  next  iterate  those 

matching  values  of 

yy'  =  y0(q)  ;  %  xO  and  yO 

WW  =  wO ; 
z  z  =  z  0  ; 
pause ; 

plot  'xx, yy ,  ' g . ' , wO , zO , ' w+ ' , xline, yline, ' w- ' ) ; hold  on; 
plot (o , s , ' w . ' , u, v, ' w . ' , e, f , ' w. ' , g, h, ' w . ' ) ;  hold  off; 
title  ( 'Movement  of  points  in  grid  with  sam.e  bina'ry 
sequence ' ) ; 
pause ; clg 

end 

axis  ; 

vec  =  ( 1 /2 )."(( 1 : length (qvec) )-l); 
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plot { ( 1 : length ( qvec ) ) -1 , qvec , ' r ' , ( 1 : length ( qvec ) ) -1 , qvec ( 1 i ' 
vec , ' b ' ) ; 

title ( 'Decrease  in  #  of  grid  pts  with  same  binary  sequence 
as  iter  incr' ) 
pause ; clg 

%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  original 
xvec  and  yvec . 

axis ( [-1.32  1.32  -.6  .5]); 

plot ( o, s , ' w . ' , u, V, ' w . ' , e, f , ' w . ' , g, h, ' w . ' ) ; hold  on; 
ploc (xvec ( indices ) , yvec (indices), 'gx',W,Z, 'r*',x,y, 'b.',xlin 
e,  :;^'ii.ne,  '  w'  )  ; 

title('W,Z  and  the  grid  point  with  the  same  binary 
sequence' ' ; pause, • 
held  off; 

xnew  =  xvec ( indices ) ; 
ynew  =  yvec ( indices ) ; 
axis(  [-1.32  1.32  -.6  .5]  )  ; 

plct(xline,yline, 'w',x,y, 'b.',o,s, 'w.',u,v, 'w.',e,f, 'w.',g,h 
,  'w.  '  )  ; 

hold  on; 
xnewt  =  xnew; 


ynewt  =  ynew; 
wt  =  w; 
z  t  =  z  ; 


'or  n  = 

1: 

:  i  -  1 

xOt 

= 

1.0  -  1.4  *xnewt  .''2  +  ynewt  ; 

yOt 

= 

. 3  * xnewt ; 

wOt 

= 

1.0--  1 . 4*wt .  ''2  +  zt  ; 

zOt 

= 

. 3  *wt ; 

wt 

= 

wOt  ; 

zt 

zOt  ; 

xnewt 

= 

xOt  ; 

ynewt 

= 

yOt ; 

plot (xOt ,yOt , ' gx' , wOt , zOt , ' r* ' ) ;hold  on; 
title { 'Movement  of  W, Z  and  matching  sequence  pt  until 
divergence' ) ; pause; 

plot (x0t,y0t,  'ix' ,wOt, zOt,  'i*'  )  ; 

end 

hold  off; 


U, 


%  function  [xant,  yaat]  >  grdcoapB (sp.M, Z) 

%  This  function  takes  a  point iW, Z i from  the  RT  quad  and 
Iterates  the  point  n 

%  times  in  order  to  ensure  the  point  is  on  the  Henon 
attractor; ensuring 

%  that  the  nth  iterate  is  in  the  left  quadrant).  It 
initiates  a  grid  field 

%  based  on  a  certain  spacing,  then  iterates  those  grid  field 
tints  that 

match  the  binary  string  of  (wnew, znew)  [use  (-l,-.25)  fcr 
,  Z .  .^s  long  as 

?  the  strings  match,  points  are  iterated  using  the  Henon 
recurrence.  Only  that 

%  I  these)  points  that  completely  match  are  finally  plotted 
according  to  the 
%  following  scheme: 

%  b.  original  grid  field  points 

%  r’  original  lW,Z) 

%  go  iterated  field  points  with  matching  binary  string 

at  nth  Iteration 

%  w-  Iterated  (W,Z)  at  nth  iteration 

%  gx  original  xvec,yvec  point (s)  which  matches  binary 

string  cf  W, Z 

%  We  also  show  how  W, Z  and  xvec ( indices ), yvec ( indices ) 

"walk"  to  the  Henon 

%  attractor  at  each  iteration  (with  the  attractor  on 
screen ) . 

c  =  linspace (-1.33,1.32,500)  ; 
s=-. 1033*0  +  .276; 
u  =  linspace(l. 32, 1.245, 500)  ; 
v=3.64*u  -4.6718; 
g  =  linspace(-1.06, 1.245, 500)  ; 
h=.1533*g  -.3344; 
e=l inspace (-1.06,-1.33,500); 
f=-3.407*e  -4.1119; 
spilt  =  .4098; 

a  =  (-.1083  -l;-3.64  1;-.1562  1]; 
b  =  2'760;-4. 6718;-. 3344]  ; 

c  =  [split: sp: 1.32]  ; 
d  -  • - . 6 : sp :  . 5 ]  ; 

[x,  y]  =  henreal ( 750 , -1 . 0 , -  .  25 )  ; 

xline  =  [split  split];  ,  yline  =  [-1.32  1.32]; 

lenc  =  ler.gth(c);  lend  length{d); 

xvec  =  zeros ( 1 , lenc*lend) ;  yvec  =  zeros ( 1 , lenc* lend) ; 
k  =  C  ; 

for  p  =  1 : lenc 

for  r  =  l;lend 
k=k+l ; 

if  max ( { a* [c ( p ) ; d { r ) ] )  <  b)  ==  0 
xvec(k)  =  c(p};  yveefk)  =  d(r); 
el;.,e 
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xvec(k)  =  -10;  yvec(k)  =  -10; 

end 

end 

end 

xvec  =  reshape ( xvec , lenc* lend, 1 ) ; 
yvec  =  reshape (yvec, lenc*lend, 1 ) ; 

m  =  find(xvec-=  -10  I  yvec~=  -10); 
xvec  =  xvec(m); 
yvec  =  yvec (m) ; 

wattr  =  w;  %  This  preserves  the  values  of  W,Z 

z  a  1 1  r  =  z  ; 

for  j  =  1:40  %  This 

ensures  water, zartr  is 

wattr{T+l)  =  1.0  -  1 . 4»wattr  ( j  )  ^2  zattr(j);  %  on  the 
attractor 

zattr(j+l)  =  .3*wattr(j); 

if  j  >  20  &  wattr(j+l)  >  split,  break; 
end 

end 


wnew 

=  wa  1 1  r  (  j  1 )  ; 

%  preserves  value  on  attracter 

wattr. 

zattr 

znew 

=  zattr  ( j  -t-1 )  ; 

wattr 

=  wnew; 

zattr 

=  znew; 

XX 

=  xvec ; 

%  preserves  the  values  of  xvec  £c  yvec 

yy 

=  yvec; 

axis  (  [ 

-1.32  1.32  -.6 

.51  )  ; 

indices  =  1 : length (xx) ;  %  [1  2  3  4  ..  .1 

1  =  1; 

qvec  =  length(xx); 
q  =  indices; 

lenq  =  length(xx)  %  initialization  of  q 

plot (xvec, yvec , 'b. ',wattr,zattr, 'r*', xline, yline , ' w- ' j  ; hold 
on  ; 

plot{o,s,'w.',u,v,'w.',e,f,'w.',g,h,'w.');  hold  off; 
title! 'Grid  field  and  point  chosen  on  attractor pause ; clg ; 
while  length (q)  >  1; 

axis!  [-1.32  1.32  -.6  .5] )  ; 

xO  =  yy  +  ones ( size (yy ) )  -  1.4*xx."2; 

yO  =  .3*xx; 

wO  =  znew  +  1  -  1 . 4*wnew. '"2  ; 
zO  =  .3*wnew; 

if  wO  >  split 

qtemp  =  findlxO  >  split  I  xO  ==  split); 
else 

qtemp  =  find(x0  <  split); 
end; 

if  length (qtemp)  ==  0, break ; 

else  q  =  qtemp; 

end 
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%  counts  Iterates  where  at  least 


1  =  i  +  1  ; 
one  point  matches 

qvec  =  [qvec, length(q) ] ; %  plot  of  how  #  with  same  binaiA' 
sequence 

lenq  =  length(q) 

indices  =  indices (q);  %  decreases  each  time  through  loop 

XX  =  xO(q);  %  preserves  to  next  iterate  those 

matching  values  of 

•/y  ^  yO  !q)  ;  %  xO  and  yO 

wnew  =  wO ; 
znew  =  zO; 
pause; 

plot  (XX,  y:^',  '  g  .  '  ,  wO  ,  zO  ,  '  w+  '  ,  xline,  yline,  '  w-  '  )  ;  hold  on; 
plot ( o , s , ' w . ' , u, V, ' w . ' , e, f , ' w. ' , g, h, ' w. ' ) ; 
tit le ( 'Movement  of  points  in  grid  with  same  binary 
sequence ' ) 

pause;  cig;hold  off; 

end 

axis([l  2  3  4]);axis; 

vec  =  ( 1 /2 ). ^ (( 1 : length (qvec) ) -1 ) ; 

plot (  ( 1 : length ( qvec )  ) -1 , qvec ,  ' r ' ,  ( 1 : length ( qvec )  ) -  1 , qvec ; 1 ■ * 
vec ,'b'); 

t  It  le  (' Decrease  in  #  of  grid  points  with  same  binar'/ 
sequence  vs  l/2"n') 
pause; clg; 

%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  original 
xvec  and  yvec . 

axis ([-1.22  1.32  - . 6  .5]); 

plot(xline,yline, 'w',x,y, 'b.',o,s, 'w.',u,v, 'w.',e,f, 'w. ',g,h 
, 'w. ' ) ;hold  on; 

plot (xvec ( indices ) , yvec (indices ) ,  ' gx ' , wat  tr , zatt r ,  ' r’ '  ; held 
on; 

title ( 'wattr, zattr  &  grid  pt  with  same  binary 
sequence' ) ; pause ; 
hold  off; 

xnew  =  xvec ( indices ) ; 

ynew  =  yvec ( indices ) ;  %  now  we  show  a  plot  of  how 
wattr, zattr  and  the  point 

%  that  generates  the  same  binary  sequence 

xnew, ynew 

%  walk  around  the  attractor 
axis( [-1.32  1.32  -.6  .5] ) ; 

plot (xline,yline, 'w' ,x,y, 'b. ' ,o,s, 'w. ' ,u,v, 'w. ' ,e, f, 'w. ' ,g,h 

,  'w.  '  )  ; 

hold  on; 

xnewt  =  xnew; 

ynewt  =  ynew; 

wattrt  =  wattr; 


for  n  =  1 : 1  -  1 

xOt  =  1.0  -  1.4*xnewt.'2  >  ynewt ; 
yOt  =  .3*xnewt; 

wOt  =  1.0  -  1 . 4*wattrt . "2  +  zatcrt; 
zOt  =  .3»wattrt; 
wattrt  =  wOt; 
zattrt  =  zOt; 
xnewt  =  xOt ; 
ynewt  =  yOt; 

plot (xOt , yOt , ' gx' , wGt , zOt , ' r’ ' ' ; 
title i 'Movement  of  pt  chosen  on  HA  and  matching  sequence 
until  diverge'); 

pause;  plot (xOt , yOt ,  ' ix' , wOt , zOt ,  ' i*  '  )  ; 
e.nd 

held  off; 


%  f\inction  [xnat,  ymat]  «  grdcoaip4  (sp,  W,  Z) 

%  ONLY  DIFF  BT/y'N  3i4  IS  THAT  4  WEEDS  QUADR,  3  USES  WHOLE 
RECTANGLE . 

%  This  function  takes  a  point  (W,  Z '  froti  the  quadr  and 
Iterates  the  point  n 

%  times  in  order  to  ensure  the  point  is  on  the  Hencn 
at  tractor ( ensuring 

^  that  the  nth  iterate  is  in  the  left  quadrant) .  It 
initiates  a  grid  field 

*  based  on  a  certain  spacing,  then  iterates  those  grid 
Clints  that 

t  match  the  binary  string  of  v’^ieW,Znew.  [use  (-I, -.25' 
W , Z : .  As  long  as 

«  the  strings  match,  points  are  iterated  using  the  Hem 
recurrence.  Only  that 

*  those  points  that  com.pletely  match  are  finally  plot 
according  to  the 

«  fill  owing  s  chem.e  : 

¥  b.  original  grid  field  points 

¥  r*  original  'WLZ) 

¥  go  Iterated  field  points  with  matching  binary'  s 

at  nth  Iteration 

*  •w+  Iterated  (W,Z)  at  nth  iteration 

¥  gx  original  xvec,yvec  point  (s'  which  m.atches  bi 

string  of  W, Z 
format  long 

¥  We  also  sho’w  ho-w  W,Z  and  xvec  ( indices  ),  yvec  ■' indices 
"walk"  to  the  Henon 

i  attractor  at  each  iteration  (with  the  attractor  on 
screen''  . 

o=lin3pace (-1.33,1.32,500); 

3=-. 1083*0  +  .276; 
u  =  linspace(1.3  2, 1.245, 50  0)  ; 
v=3.64*u  -4.6718; 
g  =  linspace (-1.0  6,1.245,500)  ; 
h=.1533*g  -.3344; 
e  =  linspace (-1.06, -1.33,500)  ; 
f=-3.4C7*e  -  4.1119; 
sclit  =  .4098; 

a  =  [3.4074  1;  -.1083  -1;  -.1562  1]; 
b  =  [-4.1119;  -.2760;  -.3344]  ; 

C  =  [ -1 . 32 : sp : split ] ; 

d  =  [-.6:sp:.5j; 

if  min((a*[w;z])  >  b)  ==  0  I  w  ^  split 

disp (' initial  value  is  not  in  left  quadr' 
return ; 

end 

[x,  y]  =  henreal (750 , -1 . 0 , - . 25 ) ; 
xline  =  [split  split];  yline  =  [-1,32  1.32]; 

lenc  =  length(c);  lend  =  length(d); 

xvec  =  zeros  ( 1 ,  lenc*lend',  ;  yvec  =  zeros  ( 1 ,  lenc* lend'  ,- 
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k  =  0; 

fcr  p  =  l:lenc 

for  r  =  l:lend 
k  =  k  +  1; 

%  This  next  section  puts  everything  outside  of  the  quad  to 

i  -10, -10) 

if  max((a  *  [c (p) ;d(r) ] )  <  b)  ==  0 
xvec(k)  =  c(p);  yvec (k) =d(r) ; 
else 

xvec(k)  =  -10;  yvec(k)  =  -10; 

end 

end 

end 

xvec  =  reshape (xvec , lenc* lend, 1 ) ; 
yvec  =  reshape (yvec , lenc* lend, 1 ) ; 

m  =  findlxvec  ~=  -10  I  yvec  ~=  -10); 
xvec  =  xvec(m); 
yvec  =  yvec (m) ; 

wattr  =  w;  %  This  preserves  the  values  of  W, z 

z  a  1 1  r  =  z  ; 

for  ]  =  1:40  % 

ensures  wattr, zattr  is 

wattr(j+l)  =  1.0  -  1 . 4*wattr  ( j ) ''2  +  zattrlj);  %  on  the 
attractor 

zattr  (j  +  1)  =  .3*wattr(:j); 

if  j  >  20  &  wattr{j+l)  <  split,  break; 
end 

end 

wnew  =  wat  t  r ( j  + 1 ) ; 
znew  =  zattr ( j+l ) ; 
wattr  =  wnew; 
zattr  =  znew; 

XX  =  xvec;  %  preserves  the  values  cf 

xvec  Sc  yvec 
yy  =  yvec; 

axis (  1-1.32  1.32  -.6  .5]  )  ; 

indices  =  1 : length (xx) ;  %  [1  2  3  4  .  .  .] 

1  =  1; 

qvec  =  length(xx); 

q  =  indices; 
lenq  =  length (xx) 

plot ( xvec , yvec , ' b . ' , wattr , zattr , ' r* ' , xline, yline, ' w' ) ; hold 
on; 

plot ( o, s , ' w . ' , u, V, ' w . ' , e, f , ' w. ' , g, h, ' w . ' ) ;  hold  off; 
title('Grid  field  and  point  chosen  on  attractor '); pause ; clg ; 
while  length (q)  >  1; 

axi3( [-1.32  1.32  -.6  .5] )  ; 

xG  =  yy  +  ones  ( size  (yy)  )  -  1.4*xx.''2; 

yO  =  .3*xx; 

wO  =  znew  +  1  -  1 , 4*wnew.  "'2 ; 
zO  =  .3*wnew; 
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•  •  •  • 


•  • 


•  • 


if  wO  >=  split 

qtemp  =  find(xO  >  split  I  xO  ==  split); 
else 

qtemp  =  findlxO  <  split); 
end; 

if  length (qtemp)  ==  0, break; 
else  q  -  qtemp; 
end  ; 

1  =  11-1;  %  counts  Iterates  where  at  least 

one  point  matches 

qvec  =  [qvec ,  length  (q)  ];  %  plot  of  how  #  with  same  binairy 
sequence 

lenq  =  length (q) 

indices  =  indices (q);  %  decreases  each  time  through  loop 

XX  =  x0(q);  %  preserves  to  next  iterate  those 

miatching  values  of 

yy  =  yO  iq)  ;  %  xO  and  yO 

wnew  =  wO ; 
cnew  =  zO; 
pause; 

plot  (XX,  '  g .  '  ,  wO  ,  zO ,  '  w-t- '  ,  xiine,yline,  '  w'  )  ; hold  on; 
plot(o,s, 'w.',u,v, 'w.',e,f, 'w.',g,h, 'w.') ; 
t itle (' Movement  of  points  in  grid  with  same  binary 
sequence ' ) ; 

pause;  clg;  held  off; 

end 

axis ( [ 1  2  3  4 ] ) ;  axis ; 

vec  =  (l/'2).''((l;length(qvec))-l); 

plot  !  ( 1 :  length  (qi/ec )  )  -1 ,  q’/ec,  '  r '  ,  ( 1 :  length  ( qvec )  )  -1 ,  qvec  ( 1 ) 
vec , ' b '  ;  ; 

title  ( 'Decrease  in  #  of  grid  points  with  same  binary' 
sequence  vs  l'2^n') 
pause ; clg ; 

%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  origina 

xvec  and  yvec . 

axis ([-1.32  1.32  -.6  .5]); 

plot ixline,yline, 'w' ,x,y, 'b. ' ,o,3, 'w. ' ,u,v, 'w. ' ,e, f, 'w. ' ,g, 

, ' w . ' ) ; hold  on ; 

plot (xvec ( indices) ,yvec (indices) , 'gx' ,wattr, zattr, ' r* ' ) ; 
%hcld  on; 

t itle i 'wattr, zattr  &  grid  pt  with  same  binary 
sequence' ) ; pause ; 
hold  off; 

xnew  =  xvec ( indices) ; 

ynew  =  yvec ( indices ) ;  %  now  we  show  a  plot  of  how 
wattr, zattr  and  the  point 

%  that  generates  the  same  binary  sequence 

xnew, ynew 

%  walk  around  the  attractor 


axis ( [ -1 

.32  1.32  -.6 

plot (x,y, 'b. ' ,xline, 

h, 'w. ' ) ; 

hold  on; 

xnewt  = 

xnew  ; 

ynewt  = 

ynew  ; 

wattrt  = 

wattr; 

zattrt  = 

zattr; 

for  n  = 

1:1-1 

xOt 

=  1.0  -  1.4* 

yOt 

=  . 3  *xnewt ; 

wOt 

=  1.0  -  1.4* 

z  0 

=  .3*wattrt; 

wattrt  = 

wOt  ; 

zattrt  = 

zOt  ; 

xnewt  = 

xOt  ; 

^*n0wr  = 

yOt  ; 

plot 

(xOt ,y0t , ' gx 

t ir ie I ' Movement  of  pt  chosen  on  HA  and  matching  sequence  pt 
until  diverge' ) ; 

pause ; plot (x0t,y0t,  'ix' ,wOt, zOt,  'i*'  ) ; 
end 

held  off; 


%  function  [xmat,  ynat]  s  proof  Kspcrs.  spfn,  incr) 

%  THIS  FUNCTION  GOES  BACK  TO  *»GRDC0MP1*’  AND  USES  IT  OVEP. 
MANY  SP  VALUES . 

%  This  function  takes  a  point (W, Z) from  the  RT  quad, init lates 
%a  quad  grid  field  based  on  a  certain  spacing,  then  iterates 
%tho3e  grid  field  points  that  match  the  binary  string  of 
% ■ W, Z ) [use ( . 1 , . 1 )  for  now.  As  long  as  the  strings  match, 
Ipcints  are  iterated  using  the  Henon  recurrence.  Only  that 
%; those ) points  that  completely  match  are  finally  considered, 
w  =  .5; 
z  =  .1; 

split  =  .4098; 

for  sp  =  spf n : incr : spcrs ,  %  MAJOR  OUTER  LOOP 

sp  =  sp 

a  =  [-.1083  -l;-3.64  1;-.1562  1 ]  ; 
b  =  [-. 2760  ;-4. 6718;-. 3344]  ; 
c  =  [split : sp : 1 . 32 ]  ; 
d  =  [ - . 6 : sp :  . 5 ]  ; 

lenc  =  iength(c);  lend  =  length(d); 

xvec  =  zeros ; 1 , lenc* lend) ;  yvec  =  zeros ( 1 , lenc’lend) ; 
if  min((a*[w;z])  >  b)  ==  0  I  w  <  split 

disp i ' init lai  value  is  not  in  the  right  quadr'i 
return ; 

end 
k  =  0  ; 

f or  p  =  1 : lenc 

for  r  =  l:lend 
k=k+l; 

if  max ( (a* [c (p) ; d ( r ) 1 '  <  b)  ==  0 
xvec{k)  =  c(p);  yv€C(k)  =  d{r); 
else 

xvec(k)  =  -10;  yvec(k)  =  -10; 

end 

end 

end 

xvec  =  reshape (xvec, lenc*lend, 1) ; 
yvec  =  reshape (yvec, lenc*lend, 1) ; 

m  :::  find(xvec  ~=  -10  I  yvec  -=  -10); 

.xvec  =  xvec(m);  yvec  =  yvec(m); 

.XX  =  xvec ;  yy  =  yvec  ; 

•/ZW  =  w ;  Z  Z  =  Z  ; 

indices  =  1 : length (xx) ;  %  [1234  .  .  .] 

1  =  1; 

qvec  =  length (xx); 

q  =  indices; 
lenq  =  length(xx); 
while  length (q)  >  1; 

xO  =  yy  +  ones (size (yy) ) -1 ,4*xx. ^2; 
yO  =  .3*xx; 

wO  =  zz  +  1  -  1.4*ww.^2; 
zO  =  .3  *ww; 
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if  wO  >  split 

qtemp  =  findfxO  >  split  I  xO  ==  split); 
else 

qtemp  =  find(xO  <  split); 
end; 

if  length (qtemp)  ==  0, break; 

else  q  =  qtemp; 

end 

i  =  1  +  1;  %  counts  iterates  where  at 

least  one  point  matches 
lenq  =  length { q) ; 

qvec  =  [qvec, length (q) ]; %  plot  of  how  #  with  same  binary 
sequence 

indices  =  indices (q);  %  decreases  each  time  through  Icop 

XX  =  xO(q);  %  preserves  to  next  iterate  those 

matching  values  of 

y'y  =  yO{q)  ;  %  xO  and  yO 

•ww  =  wO  ; 
zz  =  zO; 

end 

%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  origina 

xvec  and  yvec . 

xnew  =  xvec { indices } ; 

•/new  =  yvec  ( indices )  ;  %  now  we  show  a  plot  of  how  W,Z  and 

%the  point  that  generates  the  same 
%binary  sequence  xnew,ynew  walk  around  the  attractor 
xnewt  =  xnew; 
ynewt  =  ynew; 
wt  =  w  ; 
z  t  =  z  ; 

for  n  =  l:i  -  1 

xOt  =  1.0  -  1 . 4*xnewt . ''Z  +  ynewt; 
yOt  =  .3*xnewt; 
wOt  =  1.0  -  1.4*wt.^2  +  zt; 
zOt  =  .3*wt; 
wt  =  wOt; 
z  t  =  z  0 1  ; 
xne'wt  =  xOt; 
ynewt  =  yOt; 
end 

nvec  =  [nvec,n]; 
spvec  =  Lspvec,sp]; 
clear  xvec; clear  yvec; 

end  %  END  OF  MAJOR  LOOP 

axis  ; 

plot ( spvec , nvec ,'w*'); 

title ('SP  VS.  ITERATIONS  UNTIL  DIVERGENCE'); 


%  function  [xaat,  yaat]  «  proof 2 (spcra, apfn, Incr) 

%  THIS  FUNCTION  GOES  BACK  TO  •*GRDC0MP2**  AND  USES  IT  OVER 
MANY  SP  VALUES. 

%  This  function  takes  a  (W,Z)  from  the  left  quad, initiates  a 
quad  grid  field 

%  based  on  a  certain  spacing,  then  iterates  those  grxd  field 
points  that 

%  match  the  binary  string  of  (W,Z)  [use  (-l,-.25)  for  now]. 
.i.s  long  as  the 

%  strings  match,  points  are  iterated  using  the  Henon 
.recurrence.  Only  that 

%  i those;  points  that  completely  match  are  finally 
co.nsidered . 

w  =  ~  1  .  0  ; 
z  3  -.25; 
split  =  .4098; 

for  sp  =  spf n ; incr : spcrs ,  %  MAJOR  OUTER  LOOP 

sp  =  sp 

a  [3.40'^  1;  -.1083  -1  ;  -.1562  1]; 
b  =  [-4.1119;  -.2760;  -.3344]; 

C  =  [ -1 . 32 : sp : split ] ; 

d  -  [ - . 6 : sp ; . 5 ] ; 

lenc  =  length(c);  lend  =  length(d); 

xvec  =  zeros ( 1 , lend* lenc ) ;  yvec  =  zeros ( 1 , lenc* lend) ; 
if  min ( ( a* [w; z ] )  >  b)  ==  0  1  w  >  split 

disp  i  ' init ial  value  is  not  in  the  left  quadr ' ) 
return; 

end 

k  =  0  ; 

for  p  =  1 ; lenc 

for  r  =  l;lend 
k=k+l ; 

1 f  max ((a*[c(p);d(r)])  <  b)  ==  0 
xvec(k)  =  c(p);  yvec(k)  =  d(r); 
else 

xvec(k)  =  -10;  yvec(k)  =  -10; 

end 

end 

end 

xvec  =  reshape (xvec, lenc*lend, 1 ) ;  %makes  a  column  vector 
yvec  =  reshape (yvec, lenc*lend, 1) ; 

m  =  find ( xvec ~=  -10  I  yvec~=  -10); 

xvec  =  xvec  f.m)  ;  yvec  =  yvec'm); 

XX  =  xvec ;  yy  =  yvec ; 

ww  =  w ;  z  z  =  z  ; 

indices  =  1 : length (xx) ;  %  [1234  .  .  .] 

1  =  1; 

qvec  =  length(xx); 

q  =  indices;  %  initialization  of  q 

lenq  =  length(xx); 
while  length (q)  >  1; 
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xO  =  yy  +  ones  ( size  (yy)  )  -  1.4»xx.''2; 
yO  =  .3*xx; 

wO  =  zz  1  -  1.4*ww."2; 
zO  =  .3*ww; 

if  wO  >  split 

qtemp  =  findfxO  >  split  !  xO  ==  split); 
else 

qtemp  =  find(xO  <  split); 
end; 

if  length (qtemp)  ==  0, break; 

else  q  =  qtem.p; 

end 

1=1+1; 
lenq  =  length(q); 
qvec  =  [qvec , length (q) ] ; 
indices  =  indices (q); 

XX  =  xO(q);  %  preserves  to  next  iterate  those 

matching  values  of 

yy  =  yO (q) ;  %  xO  and  yO 

ww  =  wO ; 
z  z  =  z  0  ; 

end 

%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  original 
xvec  and  yvec . 

xnew  =  xvec ( indices ) ; 
ynew  =  yvec ( indices ) ; 
xnewt  =  xnew; 
ynewt  =  ynew; 
wt  =  w, 
z  t  =  z; 

forn=l:i-l 

xOt  =1.0  -  1.4*xnewt.^2  +  ynewt; 
yOt  =  .3*xnewt; 
wOt  =  1.0  -  1.4*wt.^2  +  zt; 
zOt  =  .3*wt; 
wt  =  wOt ; 
z  t  =  z  0 1  ; 
xnewt  =  xOt; 
ynewt  =  yOt; 
end 

nvec  =  [nvec,n]; 
spvec  =  [3pvec,sp]; 
clear  xvec; clear  yvec; 

end  %  END  TO  MAJOR  OUTER  LOOP 

plot ; spvec, nvec, 'w* ' ) ; 

title('SP  VS.  ITERATIONS  UNTIL  DIVERGENCE’); 
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»-ti  (D  U*  (T)  ft,  ij  3-‘  o>  «  )  o\‘  'iv  >]  ‘.-v  'O  -A*'  U>  •A'^  'AS'  'V  O  ^  o^P 


fTinction  [xmat,  ymat]  *  proo£3  (spcrs,  apfn,  incr) 

THIS  FXN  USES  ( W, Z ) = t - 1 , - . 2 5 . .  ONE  ATTR  PT  IS  COMPARED 
VER  MANY  SP  VALUES 

ONLY  DIFF  BTWN  3i4  IS  THAT  4  WEEDS  QUADR,  3  USES  WHOLE 
ECTANGLE . 

This  function  takes  a  point (W, Z ; from  the  quadr  and 
terates  the  point  n 

times  ir.  order  to  ensure  the  point  is  on  the  Hencn 
t tract  or  ensuring 

that  the  nth  iterate  is  in  the  left  quadrant) .  It 
nitrates  a  grid  field 

based  on  a  certain  spacing,  then  iterates  those  grid  field 
tints  that 

m.aton  the  binary  string  of  (wnew,  znew)  [use  (-l,-.25  ror 
, Z ; .  As  long  as 

the  strings  match,  points  are  iterated  using  the  Hencn 
ecurrence.  Only  that 

■tricse'  points  that  completely  match  are  finally 
'A  =  .  b  ; 

plit  =  .4C9S; 

attr  =  w;  %  This  preserves  the  values  of  W,  Z 

a  1 1  r  =  z  ; 

or  :  =  1:40  « 

nsures  wattr,zattr  is 

wattr::-^!)  =1.'^  -  1 . 4*wattr  ( j  )  ^2  +  zattr(j);  %  on  the 
ttractcr 

zattr't^l)  =  . 3 *wattr ( j ) ; 

if  :  :•  20  St  wattrij  +  l)  >  split,  break; 

end 
nd 

new  =  wa  1 1  r  ^  j  1 )  ; 
new  =  zattr ( 3  +1 )  ; 

or  sp  =  spf n : incr ; spcrs ,  %  MAJOR  OUTER  LOOP 

sp  =  sp 
wattr  =  wnew; 
zattr  =  znew; 

a  =  :-.10b3  -1;  -3.64  1;  -.1562  1]; 
b  =  [- .2760  ; -4.6718; - .3344]  ; 

C  =  [ spilt : sp ; 1 . 32 ]  ; 

d  =  [ - . 6 : sp :  .  5 ]  ; 

lenc  =  length(c);  lend  =  length(d); 

.xvec  =  zeros  ( 1 ,  lenc*lend)  ;  yvec  =  zeros  { 1 ,  lenc’lend :  ; 
if  .min  (  ( a*  [w;  z ]  )  >  b/  ==  0  i  w  <•  split 

disp (' initial  value  is  not  in  right  quadr'' 
return ; 

end 

k  =  C  ; 

for  p  =  1 : lenc 

for  r  =  1 : lend 
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k=k+l ; 

%  This  next  section  puts  everything  outside  of  the  quad 
to  (-10,-10) 

if  max ((a  •  [c(p);d(r)])  <  b)  ==  0 
xvec(k)  =  c{p);  yvecik}  =  d(r); 
else 

xvec(k)  =  -10;  yvec(k)  -  -1C; 
end 

end 

xvec  =  reshape (xvec, lenc* lend, 1 ) ; 
yvec  -  reshape (yvec, lenc*lend, 1 ) ; 

m  =  find (xvec  ~=  -10  i  yvec  -=  -10); 

xvec  =  xvec{m);  yvec  -  yvec(m); 

xx  =  xvec ;  vy’  =  yvec  ; 

indices  =  1 : length (xx) ;  %  [1234  .  .  .] 

1  =  1; 

qvec  =  length(xx); 
q  =  indices; 
lenq  =  length(xx); 
while  length (q)  >  1 ; 

xO  =  yy  +  ones ( size (yy) )  -  1.4*xx."2; 
yO  =  .3*xx; 

wO  =  znew  1  -  1.4*wnew.^2; 
zO  =  .3*wnew; 

if  wO  >=  split 

qtemp  =  find(x0  >  split  I  xO  ==  split); 
else 

qtemp  =  find(x0  <  split); 
end; 

if  length (qtemp)  ==  0, break; 

else  q  =  qtemp,- 

end; 

i  =  1  1 ;  %  counts  iterates  where  at  least 

one  point  matches 

qvec  =  [qvec,  length  (q)  1  ;  %  plot  of  how  *  with  sam.e 

binary  sequence 

lenq  =  length (q); 

indices  =  indices (q);  %  decreases  each  time  through  loop 

x:-:  -  xC(q);  %  precor”es  fo  nex*"  iterate  those 

matching  values  of 

yy  =  yO (q) ;  %  xO  and  yO 

wne'w  =  wO  ; 
znew  =  zO; 

end 

%  In  this  way,  we  need  the  index  key  'indices'  '-c  '■el''  to 
which  grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  original 

xvec  and  yvec . 

xnew  =  xvec (indices); 


♦ 


•  • 


•  •  • 
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ynew  =  yvec ( indi ces ) ;  %  now  we  show  a  plot  of  how 
wattr.zattr  and  the  point 

%  that  generates  the  same  binary  segue 

xnew , ynew 

%  walk  around  the  attractor 


xnewt  =  xnew; 
ynewt  =  ynew; 
wattrt  -  wattr; 
cattrt  =  cattr; 
f ;  r  n  =  1  :  i  -  1 

xOt  =  1.0  -  1.4 'xnewt. '2  *  ynew:  ; 
yOt  =  . 3 'xnewt; 

w't  =  1.0  -  1 . 4'wattrt  .  "2  cattrt; 
c  ■:t  .3  'wattrt  ; 

wattrt  =  wOt ; 
cattrt  =  cOt; 


.xnewt  =  xOt 


ncec  =  [nvec,n’ ; 
spvec  =  [3pvec,sp]; 
clear  xvec;  clear  yvec ; 

end  %  END 


>!AJ 


axis  ; 

plot ' spvec, nvec, 'w* ' ) ; 

title' 'S?  VS.  ITERATIONS  UNTIL  DIVERGENCE'’; 


ynew  =  yvec ( indices ) ;  %  now  we  show  a  plot  of  howwattr, 

%  zattr  and  the  point  that  generates 
%the  same  binary  sequence  xnew,ynew 
%  walk  around  the  attractor 

xnewt  =  xnew; 
ynewt  =  ynew; 
wattrt  =  wattr; 
cattrt  =  zattr; 
for  n  =  1:1  -  1 

xGt  =  1.0  -  1.4*xnewt."2  +  ynewt; 
yOt  =  .3*xnewt; 

wOt  =  1.0  -  1 . 4*wattrt . ''2  -t-  zattrt; 
zOt  =  .3*wattrt; 
wattrt  =  wOt; 
zattrt  =  zOt; 


xnewt  =  xOt; 
ynewt  =  yOt; 
end 

nvec  =  [nvec , n] ; 
spvec  =  [spvec,sp]; 
clear  xvec;  clear  yvec; 

end  %  E! 

axis  ; 

plot ( spvec , nvec , ' w* ' ) ; 
title('SP  VS.  ITERATIONS  UNTIL  DIVERGENCE 


%  END  TO  MAJOR  LOOP 
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%  fxmction  [xatat,  yniat]  >  prf  1  (spcrs,*pfn,  incr) 

%  This  function  cakes  thirty-one  points  (W,Z)  from  the  P.T 
quad, init lates  a 

%  quad  grid  field  based  on  a  certain  spacing,  then  icerares 
hose  grid  field 

points  that  match  the  binary  string  of  (W, 2 ) [use ( . 1 , . 1 
or  now.  As  long  as 

the  strings  match,  points  are  iterated  using  the  Henon 
ecurrence.  Only  that 

i those)  points  that  completely  match  are  finally 
onsidered . 

11  =  0 ; 

or  w  3  . 5 :  . 02  :  1  .  1  %  MAJOR  OUTER  LOCT 

c  =  .1; 

11  =  11  +  1 
plit  =  .4098; 

cr  sp  =  spf n : incr : spcrs ,  %  MINOR  OUTER  LOOP 

sp  =  sp 

a  =  [-.1083  -l;-3.64  1;-.1562  1]; 
b  =  [-. 2760;-4. 6718;-. 3344] ; 

C  =  [ split ; sp : 1  .  32 1  ; 
d  =  [ - . 6 : sp :  . 5 ]  ; 

lenc  =  length(c);  lend  =  lengthid  ; 

.xvec  =  zeros  ( 1 ,  lenc'lend!  ;  yvec  =  zeros  ( 1 ,  lenc’lend,  ; 
if  min ( ( a* [w; z ] )  >  b)  ==  0  I  w  <  split 

disp ( ' init ial  value  is  not  in  the  right  quadr ' ) 
return ; 

end 
k  =  0  ; 

for  p  =  lilsnc 

for  r  =  1 ; lend 
k=k*i ; 

if  max! (a* [c (p) ;d{r) ] )  <  b)  ==  0 
xvec(k)  =  c(p);  yvec(k)  =  d(ri; 

0  —  3  © 

xvecik)  =  -10;  yvec'k)  =  -10; 

end 

end 

xvec  =  reshape (xvec , lenc* lend,  1 )  ; 
yvec  =  reshape (yvec, lenc*lend, 1 ) ; 

m,  =  find  (xvec  ~=  -10  I  yvec  -=  -10); 

xvec  =  xvec{m) ;  yvec  =  yvec(m) ; 

xx  =  xvec ;  yy  =  yvec ; 

ww  =  w ;  z  z  =  z ; 

indices 
.  1 
qvec 
q 

lenq 


=  1 : length (xx) ;  %  [1  2  3  4  .  .  .] 

_  T  . 

=  length(xx); 

=  indices; 

=  length (xx) ; 


▼ 


♦ 


4r 
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while  length (q)  >  1; 

xO  =  yy  +  ones  ( size  (yy )  )  -  1.4*xx.''2; 
yO  =  .3*xx; 

wO  =  zz  +  1  -  1.4*ww."2; 
zO  =  .3*ww; 

if  wO  >  split 

qtemp  =  find(x9  >  split  I  xO  ==  split); 
else 

qtemp  =  findlxO  <  split); 
end; 

if  length ( qtemp )  =-  0, break; 

else  q  =  qtemp; 

end 

i  =  i  +  1;  %  counts  iterates  where  at 

least  one  point  matches 
lenq  =  length ( q) ; 

qvec  =  [qv'ec,  length  (q)  ];  %  plot  of  how  #  with  same  binary' 
sequence 

indices  =  indices (q);  %  decreases  each  time  through  loop 

XX  =  xO(q);  %  preserves  to  next  iterate  those 

m.atching  values  of 

-  y0{q) ;  %  xO  and  yO 

ww  =  wO ; 
z  z  =  z  0  ; 

end 

%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  original 

xvec  and  yvec . 

xnew  =  xvec ( indices) ; 

ynew  =  yvec ( indices ) ;  %  now  we  show  a  plot  of  how  W, Z  and 
the  point 

%  that  generates  the  same  binary  sequence 

xnew, ynew 

%  walk  around  the  attractor 


xnewt  =  xnew; 
ynewt  =  ynew; 
wt  =  w  ; 
z  t  =  z; 

forn=l:i-l 

xOt  =  1.0  -  1.4*xnewt.^2  +  ynewt; 
yOt  =  .3*xnewt; 
wOt  =  1.0  -  1.4*wt.^2  +  zt; 
zOt  =  .3*wt; 
wt  =  wOt; 
zt  =  zOt; 


xnewt 

ynewt 

end 

nvecl  = 
spvecl  = 


=  xOt; 

=  yOt; 

[nvecl , n] ; 
[spvecl, sp] ; 


•  •  • 


1 •  . . ft 
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%  END  MINOR  LOOP 


clear  xvec; clear  yvec; 
end 

nvec2  =  [nvec2  ,  nvecl  ]  ;  • 

spvec2  =  [spvec2 , spvecl ] ; 
clear  xvec;  clear  yvec; 

end  %  ENT)  MAJOR  LOOP 

axis{[0  max(spvec2)  0  max (nvec2 ) ] ) ; 
plot ( 3pvec2 , nvec2 

t;icie('SP  VS.  ITERATIONS  UNTIL  DIVERGENCE');  • 
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%  function  [xaat,  ymat]  «  prf2 (cpcrs, >pfn, incr) 

%  THIS  FUNCTION  GOES  BACK  TO  **GRDCOMP2**  ANL'  USES  IT  OVER 
MANY  SP  VALUES . 

%  (AND  OVER  THIRTY-ONE  W, Z  VALUES) 

%  This  function  takes  a  (W,Z)  from  the  left  quad, initiates  a 
quad  grid  field 

%  based  on  a  certain  spacing,  then  iterates  those  grid  field 
points  that 

%  match  the  binary  string  of  (W,Z)  [use  (-l,-.25)  for  now]. 
As  long  as  the 

%  strings  match,  points  are  iterated  using  the  Henon 
recurrence.  Only  that 

%  -those)  points  that  completely  match  are  finally 
considered . 

11  =  0  ; 

for  w  =  -  .  80  :  .  04  :  .  4  %  MAJOR  -OUTER  L-OOP 

z  =  -.01 
li  =  ii  +  1 
split  -  .4098; 

for  sp  =  spf n : incr : spcrs ,  %  MINOR  OUTER  LOOP 

sp  =  sp 

a  =  [3.407  1;  -.1083  -1;  -.1562  1]; 
b  =  [-4.1119;  -.2760;  -.3344] ; 
c  =  [ -1 . 32 : sp : split ] ; 
d  =  [ - . 6 : sp : . 5  ]  ; 

lenc  =  length(c);  lend  =  length(d); 

xvec  =  zeros ( 1 , lend*lenc ) ;  yvec  =  zeros ( 1 , lenc*lend) ; 

if  min ( ( a* [w; z ] )  >  b)  ==  0  I  w  >  split 

disp (' initial  value  is  not  in  the  left  quadr ' ) 
return ; 

end 

k  =  0; 

for  p  =  Irlenc 

for  r  =  1 : lend 
k=k+l ; 

if  max( (a* [c (p) ;d{r) I )  <  b)  ==  0 
xvec(k)  =  c(p);  yvec(k)  =  d(r); 
else 

xvec(k)  =  -10;  yvec{k)  =  -10; 

end 

end 

and 

xvec  =  reshape (xvec , lenc* lend, 1 ) ;  %makes  a  column  vector 
yvec  =  reshape (yvec, lenc*lend, 1) ; 

m  =  find(xvec~=  -10  I  yvec-=  -10); 

xvec  =  xvec(m);  yvec  =  yvec(m); 

XX  =  xvec ;  yy  =  yvec ; 

ww  =  w ;  z  z  =  z ; 

indices  =  1 ; length (xx) ;  %  [1  2  3  4  .  .  .] 

i  =  1 ; 

qvec  =  length(xx); 


q  =  indices; 
initialization  of  q 
lenq  =  length(xx); 
while  length(q)  >  1; 


xO 

=  yy 

+  ones 

(sizeiyy))  -  1.4*xx.''2; 

yO 

=  .3*xx; 

wO 

=  zz 

+  1  - 

1 . 4*ww . "2 ; 

zO 

=  .3’ 

'ww; 

if  wO  >  split 

qtemp  =  find(xO  >  split  I  xO  ==  split) ; 
else 

qtemp  =  find(xO  <  split); 
end  ; 

if  length ( qtemp )  ==  0, break ; 

else  q  =  qtemp; 

end 

i  =  i  +  1  ; 
lenq  =  length (q); 
qvec  =  [qvec, length ( q) ] ; 
indices  =  indices (q); 

XX  =  x0(q);  %  preserves  to  next  iterate  those 

matching  values  of 

yy  =  yO{q) ;  %  xO  and  yO 

ww  =  wO ; 
z  z  =  z  0  ; 

end 

%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  original 
xvec  and  yvec . 

xnew  =  xvec ( indices ) ; 
ynew  =  yvec ( indices ) ; 
xnewt  =  xnew; 
ynewt  =  ynew; 
wt  =  w; 
z  t  =  z  ; 

for  n  =  1 : i  -  1 

xOt  =1.0  -  1.4*xnewt.^2  +  ynewt; 
yOt  =  .3*xnewt; 
wOt  =  1.0  -  1.4*wt.^2  +  zt; 
zOt  =  .3*wt; 
wt  =  wOt; 
zt  =  zOt; 
xnewt  =  xOt; 
ynewt  =  yOt; 
end 


•  • 


nvecl  =  [nvecl,n]; 
spvecl  =  [spvecl,sp]; 
clear  xvec /clear  yvec; 
end 


%  END  TO  MINOR  OUTER  LOOP 
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nvec2  =  (nvec2 , nvecl ] ; 
spvec2  =  [ spvec2 , spvecl ] ; 

end  %  END  TO  MAJOR  OUTER  LOOP 

axis([0  max(spvec2)  0  max (nvec2 ) ] ) ; 
plot ( spvec2 , nvec2 

titlei'SP  VS.  ITERATIONS  UNTIL  DIVERGENCE ' ) ; 


%  function  [xnat,  yaat]  «  prf 3 (apcrs, spfn, incr) 

%  THIS  FXN  USES  THIRTY-ONE  W,  Z  VALUES.  ONE  -ATTR  PT  IS 
COMPARED  OVER  MANY  SP  VALUES 

%  ONLY  DIFF  BTWN  3&4  IS  THAT  4  USES  LEFT  QUAE^R  ,  3  USES 
RIGHT  QUADR. 

%  This  function  takes  a  points  (W,Z)  from  the  quadr  and 
Iterates  the  point  n 

%  times  in  order  to  ensure  the  point  is  on  the  Henon 
at  t  ractor ( ensuring 

%  that  the  nth  iterate  is  in  the  left  quadrant).  It 
initiates  a  grid  field 

%  based  on  a  certain  spacing,  then  iterates  those  grid  field 
points  that 

%  match  the  binary  string  of  (-^ew.znew)  [use  (-l,-.25)  for 
X,  Z ]  .  .As  long  as 

%  the  strings  match,  points  are  iterated  using  the  Henon 
recurrence.  Only  chat 

%  (those)  points  that  completely  match  are  finally 
considered . 


for  w  =  .5; .02:1.1 
z  =  .1; 

11  =  ii  +  1 


split  =  .4098; 
water  =  w ; 
zattr  =  z ; 
for  j  =  1:40 
water ( j  +1 ) 
attractor 

zattr  ;  j  +1  )■ 
if  3  > 

end 


%  This  preserves 

%  ensures  water, 
=  1.0  -  1 . 4*wattr ( 3 ) ^2  + 

=  . 3 'water ( j ) ; 

20  &  wattr(j+l)  >  split. 


%  MAJOR  OUTER  LOOP 


the  values  of  W, Z 
zattr  IS 

zattr ( j ) ;  %  on  Che 


break; 


end 

wnew  =  water (j+1); 
znew  =  zattr  { j  i-l )  ; 

for  sp  =  spf n : incr : spers ,  %  MINOR  OUTER  LOOP 

sp  =  sp 
water  =  wnew; 
zattr  =  znew; 

a  =  [-.1083  -1;  -3.64  1;  -.1562  1]; 
b  =  [-.2760;-4. 6718;-. 3344]  ; 

C  =  [split : sp : 1 . 32 ] ; 
d  =  [ - . 6 : sp :  . 5 ]  ; 

lenc  =  length(c);  lend  =  length(d); 

xvec  =  zeros ( 1 , lenc'lend) ;  yvec  =  zeros { 1 , lenc'lend) ; 
if  min { (a* [w; z] )  >  b)  ==  0  I  w  <  split 

disp (' initial  value  is  not  in  right  quadr') 
return; 

end 


k  =  0; 

for  p  =  1 : lenc 
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for  r  =  l:lend 
k=k-t-l  ; 

%  This  next  section  puts  everything  outside  of  the  quad 
to  (-10,-10) 

if  max ((a  *  [c(p);d(r)])  <  b)  ==  0 
xvec(k)  =  c(p);  yvec(k)  =  d(r); 
else 

xvec(k)  =  -10;  yvec(k)  =  -10; 
end 

end 

end 

xvec  =  reshape (xvec, lenc* lend, 1) ; 
yvec  =  reshape (yvec, lenc'lend, 1 ) ; 

m  =  find (xvec  ~=  -10  I  yvec  ~=  -10); 

xvec  =  xvec(m);  yvec  =  yvec(m); 

XX  =  xvec ;  yy  =  yvec ; 

indices  =  1 :  length  (xx)  ,-  %  (1  2  3  4  .  .  .] 

1  =  1; 

qvec  =  length (xx); 
q  =  indices; 
lenq  =  length(xx); 
while  length (q)  >  1; 

xO  =  yy  -*■  ones  ( size  (yy)  )  -  1.4*xx.''2; 
yO  =  .3*xx; 

wO  =  znew  +  1  -  1.4*wnew.^2; 
zO  =  .3*wnew; 

if  wO  >=  split 

qtemp  =  find(x0  >  split  1  xO  ==  split); 
else 

qtemp  =  find(x0  <  split); 
end; 

if  length ( qtemp )  ==  0, break; 

else  q  =  qtemp; 

end; 

i  =  i  +  1;  %  counts  iterates  where  at  least 

one  point  matches 

qvec  =  [qvec, length (q) ] ;  %  plot  of  how  #  with  same 

binary  sequence 

lenq  =  length(q); 

indices  =  indices (q);  %  decreases  each  time  through  loop 

XX  =  x0(q);  %  preserves  to  next  iterate  those 

matching  values  of 

yy  =  y0{q) ;  %  xO  and  yO 

wnew  =  wO ; 
znew  -  zO; 

end 

%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  .grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  original 

xvec  and  yvec . 

xnew  =  xvec (indices); 


1  • t 
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ynew  =  yvec ( indices ) ;  %  now  we  show  a  plot  of  how 
wattr.zattr  and  the  point 

%  that  generates  the  same  binary  sequence 

xnew, ynew 

%  walk  around  the  attractor 

xnewt  =  xnew; 
ynewt  =  ynew; 
wattrt  =  ..-attr; 
zattrt  =  zattr; 
forn=l:i-l 

xOt  =  1.0  -  1 . 4 *xnewt . '2  +  ynewt; 
yOt  =  .3*xnewt; 

wOt  =  1.0  -  1 . 4*wattrt . "2  +  zattrt; 
zOt  =  .3*wattrt; 
wattrt  =  wOt; 
zattrt  =  zOt; 
xnewt  =  xOt; 
ynewt  =  yOt; 

end  %  END  TO  MINOR  LOOP 

nvecl  =  [nvecl,n]; 
spvecl  =  [spvecl.sp]; 
clear  xvec;  clear  yvec; 

end  %  END  TO  MAJOR  LOOP 

nvec2  =  [nvec2 , nvecl ] ; 
spvec2  =  [ spvec2 , spvecl ] ; 
end 

axis ; [ 0  max ( spvecO )  0  max (nvec2 ) ] ) ; 
plot ( spvec2 , nvec2 

titlel'SP  VS.  ITERATIONS  UNTIL  DIVERGENCE'); 


%  function  [xaint,  yaat]  *  prf4(cpcrs,8pfn,  Incr) 

%  THIS  FXN  USES  THIRTY-ONE  W, Z  VALUES.  ONE  ATTR  FT  IS 
COMPARED  OVER  MANY  SP  VALUES 

%  ONLY  DIFF  BTWN  3&4  IS  THAT  4  USES  THE  LEFT  QUADR,  3  USES 
THE  RIGHT  QUADR. 

%  This  funcCion  takes  thirty  points  Z '  from  the  quad  and 
Iterates  the  point  n 

?  times  in  order  to  ensure  the  point:  is  on  the  Henon 
attracioriensuring 

%  t.hat  the  nth  iterate  is  in  the  left  quadrant)  .  It 
i.niiiates  a  grid  field 

%  based  on  a  certain  spacing,  then  iterates  those  grid  field 
points  t.hat 

%  match  the  binary  stri.ng  of  (wnew,  znew)  [use  (-l,-.25)  for 
A  ,  Z  ]  .  .2^3  long  as 

%  the  strings  match,  points  are  iterated  using  the  Henon 
recurrence.  Only  that 

%  : those)  points  that  completely  match  are  finally 
considered . 


1.  — 

f cr  w  = 

-  ’.80  :  . 

z  = 

-  .  01 

IT  — 

11  +  1 

spilt  = 

.4098; 

wa  1 1  r  = 

w; 

zattr  = 

z  ; 

f or  ]  = 

1:40 

wat 

tr  (  j  +  1 ) 

zat 

tr  (]+l) 

if  ]  > 

end 

end 

wnew  = 

wattr ( 

znew  = 

zattr ( 

for  sp 

=  spfn: 

sp 

=  sp 

wdt  t  r  = 

wnew; 

zattr  = 

znew; 

%  MAJOR  OUTER  LOOP 


%  This  preserves  the  values  of  W,Z 


%  ensures 
%wattr , zatt r  is 


%attractor 


=  . 3 'wattr ( j ) ; 


split,  break; 


%  MINOR  OUTER  LOOP 


a  =  [3.4074  1;  -.1083  -1;  -.1562  1]; 
b  =  [-4.1119;  -.2760;  -.3344]; 

C  =  [ -1 . 32 ; sp : split ] ; 

d  =  [ - . 6 : sp : . 5 ] ; 

lenc  =  length(c);  lend  =  length(d); 

xvec  =  zeros ( 1 , lenc* lend ) ;  yvec  =  zeros ( 1 , lenc* lend ) ; 
if  min ( (a* [w; z] )  >  b)  ==  0  I  w  >  split 

disp (' initial  value  is  not  in  left  quadr') 
return ; 

end 


* 


•) 


♦ 
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•  •  •  • 


•  • 


• 


•  • 


for  p  =  1 : lenc 

for  r  =  l:lend 
k=k+l ; 

%  This  next  section  puts  everything  outside  cf  the  quad 

to  (-10,-10) 

1 f  max ({a  *  [c(p);d(r)])  <  b)  ==  0 
xvec(k)  =  c(p);  yvec(k)  =  dir); 
else 

xvec(k)  =  -10;  yvec(k)  =  -10; 
end 

end 

end 

xvec  =  reshape (xvec , lenc* lend, i ) ; 
yvec  =  reshape (yvec , lenc* lend,  1 ) ; 

m  =  findixvec  -=  -10  1  yvec  -10); 

xvec  =  xvec(m) ;  yvec  =  yvec(m) ; 

xx  =  xvec ;  yy  =  yvec ; 

indices  =  1 ; length (xx) ;  %  [1234  .  .  .] 

1  =  1; 

qvec  =  length(xx); 

q  =  indices; 
lenq  =  length(xx); 
while  length (q)  >  1; 

xO  =  yy  +  ones ( size (yy) )  -  1.4*xx.^2; 
y  0  =  .  3  *  XX  ; 

wO  =  znew  +  1  -  1 .4*wnew. ''2  ; 
zO  =  .3*wnew; 

if  wO  >=  split 

qtemp  =  find(x0  >  split  I  xO  ==  split); 
else 

qtemp  =  find(x0  <  split); 
end; 

if  length (qtemp)  ==  0, break; 

else  q  =  qtemp; 

end; 

1  =  i  +  1 ;  %  counts  iterates  where  at  least 

one  point  matches 

qvec  =  [qvec, length (q) ] ;  %  plot  of  how  #  with  same 

binary  sequence 

lenq  =  length (q); 

indices  =  indices (q);  %  decreases  each  time  through  loop 

XX  =  x0(q);  %  preserves  to  next  iterate  those 

matching  values  of 

yy  =  yO (q) ;  %  xO  and  yO 

wnew  =  wO ; 
znew  =  zO; 

end 

%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  original 
xvec  and  yvec . 
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xnew  =  xvec ( indices ) ; 

ynew  =  yvec i indices ) ;  %  now  we  show  a  plot  of  how 
water, zactr  and  the  point 

%  that  generates  the  same  binary  sequence 

xilfcW ,  y  P*0w 

%  walk  around  the  attractor 

xnewt  =  xnew; 

’vT.ewt  =  ynew; 
wattrt  =  wattr; 
cattrt  =  zattr; 
for  n  =  i:i  -  1 

xCt  =  1.0  -  1.4*xnewt."2  +  ynewt ; 
yOt  =  . 3 'xnewt; 

wOr  =  1.0  -  1 . 4*wattrt . "2  +  zattrt; 
zCt  =  . 3 'wattrt; 


wa  1 1  r  t 

=  wOt 

zattrt 

=  z  0 1 

xnewi 

=  xOt 

ynewt 

=  y  Ot 

end 

nvecl  =  [nvecl.n]; 
spvecl  =  [spvecl.sp]; 
clear  xvec;  clear  yvec ; 

end  %  EMD  TO  MINOR  LOOP 

nvec2  =  [nvec2 , nvecl ] ; 
spvec2  =  [ spvec2 , spvecl ]  ; 

end  %  END  TO  MAJOR  LOOP 

axis;[0  max(spvec2)  0  max (nvec2 ) ] ) ; 
plot ( spvec2 , nvec2 

title ('SF  VS.  ITERATIONS  UNTIL  DIVERGENCE'); 


%  ONTO.M 

%  THIS  PROGRAM  REQUIRES  THE  PROGRAM  HENON  TO  RUN  PROPERLY 
%  This  program  is  a  derivative  of  a  program  used  originally 
%  in  reference  13  called  "runs.m“. 


%  This  program  is  run  using 

the  word-"oncc'‘ 

*  Tc  run 

use  "X  =  henon  '  ,  0  , 

0 , ' “ ;  then  use  "onto 

%  fer  the 

appropriate  n-tupl 

e  sought . 

X  =  hen  c  n 

;5C000G, 0, 0) ; 

n  =  4 

X 

y 

r  de  t  c  2 

=  zeros (4,1); 

Idecvez3 

=  zeros ; 8  ,  1 )  ; 

decvec4  = 

zeros (16,1)  ; 

IdecvecS 

=  zeros (32,1:  ; 

Idec'/eco 

=  zeros (64,1:  ; 

1  deevee’’ 

=  zeros  ‘,128,1;  ; 

idec'/ecS 

=  zeros (256,1)  ; 

?decvecl  =  zeros ( 512 , 1 ) ; 
^decveclC  =  zeros ; 1024 , 1 ) ; 
idecvecll  =  zercs(2048,l); 
%decvecl2  =  zeros ( 4096 , 1 ) ; 
Idecvecli  =  zeros  !  3192 ,  1 )  ; 
%decvecl4  =  zeros ' 16384 , 1 )  ; 
*decvecl5  =  zeros ( 3  27  68 ,  1 )  ; 
%decvecl6  =  zeros ( 65536 , 1 )  ; 
*decvecl7  =  zeros (131072,1); 


while  m.in  ^  decvec4  )  ==  0 
1  =  1  >  1; 

%decimal  =  [2  1]  *  x(i:i*l); 

?decvec2  decimal-^1!  =  decvec2  (decimal  +  1 )  +  1; 
Idee imal  =  [4  2  ll^x?!:!-*-!); 

%decvec3 i decimal+1 )  =  decvec3 ( decimal+1 )  +  1; 
decimal=[8  4  2  l]*x(i:i+3); 

decvec4  ( decim,al  +  l )  =  decvec4  (decimal-1 )  +  1  ; 

ldecimal=[16  842  1 ] *x ( i ; i+4 ) ; 

%decvec5 ( decimal+1 )  =  decvec5 ( decimal+1 )  +  1 ; 

%decimal= [ 32  16  8  4  2  l]*x(i:i+5); 

%decvec6 ( decimal+1 )  =  decvec6 (decimal+1 )  +  1 ; 

%decimal=[64  32  16  8  4  2  1 ] *x ( i : i+6 ) ; 

%decvec7 ( decimal+1 )  =  decvec7 (decimal+1 )  +  1 ; 

%decimal= [128  64  32  16  8  4  2  1 ] *x ( i : i+7 ) ; 
IdecvecS ( decimal+1 )  =  decvec8 (decimal+l )  +  1; 
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%decimal= [256  128  64  32  16  8  4  2  1  ]*x(i:i+8); 

%decvec9 (decimal+1 )  =  decvec9 ( decimal +1 )  +  1; 

%decimal= [512  256  128  64  32  16  8  4  2  1 ) *x ( i ; 1+ 9 ) ; 

%decveclO { decimal+1 )  =  decveclO (decimal+1 )  +  1; 

%decimal= [1024  512  256  128  64  32  16  8  4  2  l]*x(i;i+10); 
%decvecl 1 { decimal+l )  =  decvecll i decimal+1 )  +  1; 

%decimal= [2048  1024  512  256  128  64  32  16  8  4  2  1 ] ‘x i i : i + 1 1 ) ; 
%decvecl2 ( decimal+1 )  =  decvecl2 (decimal+l )  +  1; 

%decimal= [4096  2048  1024  512  256  128  64  32  16  8  4  2 

1  ]  ’  x  .  1  :  1  + 1 2  )  ; 

%decvecl3 ( decimal+1 )  =  decvecl3 (decimal+l )  +  1 ; 

%v=:8192  4096  2048  1024  512  256  128  64  32  16  8  4  2  1]; 
%decimal=v*x ( i : i  +  13 )  ; 

%decvecl4 ( decimal  +  l )  =  decvecl4 (decimal  +  l )  +  1  ; 

%v-[16384  8192  4096  2048  1024  512  256  128  64  32  16  8  4  2  Ij; 
%decimal=  v  *  x{i;i+14); 

%decvecl5 ( decimal+l )  =  decvecl5 ( decimal+l )  +  1 ; 

%v=[32768  16384  8192  4096  2048  1024  512  256  128  64  32  16  8  4 

2  1  j  ; 

%decimal  =  v  *  x(i:i+15); 

%decvecl6 ( decimal  +  1 )  =  decvecl6 (decimal  +  l )  +  1  ; 

%v=[65536  32768  16384  8192  4096  2048  1024  512  256  128  64]; 
%decimal=[v  32  16  8  4  2  1]  *  x(i:i+16); 

%decvecl7 ( decimal+l )  =  decveclV (decimal+1 )  +  1; 

end 

%  i  tells  the  sequence  length  at  which  all  the  n-tuples  were 
found . 

n  %  n-tuple  looking  for 
i  %  length  of  decvec  where  found 

length (x)  %  total  length  of  original  vector  from  Henon 


%function  xx  >  flndllO(n,xO,yO) 

%  This  program  requires  that  plot (xllO , yllO , ' . ' )  be  used 
%after  the  run.  This  program  takes  n  points  on  the  attractor 
%and  finds  those  that  give  the  sequence  110.  It  then 
%iterates  the  set  3  times.  Merely  set  the  parameters  below 


%and 

type  " 

f indllO " 

n  = 

5000  ; 

xO  = 

0; 

yO  = 

0; 

%  inputs  : 

J. 

n  = 

length 

of  desired  sequence 

% 

xO  = 

initial 

X 

% 

yO  = 

initial 

y 

%outputs : 

%  XX  =  n  by  1  binary  vector 

hold  on; 

X i 1 )  =  xO ; 
yil)  =  yO; 


0=1 inspace (-1.33,1.32,500); 
s  =  -. 1083*0  -h  .276; 
u=linspace (1.32,1.245,500)  ; 
v=3 .64*u  -  4.6~18; 
g  =  l inspace (-1.06,1.245,500)  ; 
h=.1533»g  -  .3344; 
e=linspace (-1.06,-1.33,500) ; 
f=-3.407*e  -  4.1119; 


split  =  .4098; 

axis( [-1.32  1.32  -.6  .5] )  ; 

xline  =  [split  split];  yline  =  [-1.32  1.32]; 
plot (xline, yline, 'w- ' ) ;hold  on; 
plot(o,s, 'w.',u,v, 'w. ',e,f, 'w.',g,h, 'w.') ; 
title (' Exclusion  Area  of  Points  which  Never  Give  Sequence 
1100  '  ) 
pause ; 

x(2:n)  =  zeros (n-1 , 1 ) ;  %  vectors  are  preallocated 

y(2:n)  =  zeros (n-1 , 1 ) ; 

%recursive  generation  of  points 
for  i  =  l:n; 

x(i  +  l)  =  y(i)  +  1  -  1.4*x(i)''2; 
y(i+l)  =  .3  *  x(i); 

%  convert  to  binary 

if  x(i)  <=  split 

XX  (  i  )  =  0  ; 

else 
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i 


end 


XX  ( i ) 


1; 


end 


%  XX 

decvec3  =  zeros (8,1); 
i  =  0; 


while  min ( decvec3 ) ==0 
i  =  i  +  1  ; 

decimal  =[42  1 ] *xx ( i : i+2 ) ' 
if  decimal==6,  xllO 

yllO 


end 


end 


[xllO , X ( i ) ] ; 
[yllO ,y { 1 ) j ; 


%  plot (xllO , yllO  ) 


'4r 


♦ 


•  i 
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%  SUBOUAD 

%  This  program  shows  how  the  four  sides  of  the 
subquadrilateral 

%  shift  as  the  sides  are  iterated  3  times.  This  ensures  t.hat 
any  pt . 

%  falling  in  the  subquad  will  not  give  the  sequence  HOC. 

%  To  run  type  "subquad" 

XX  =  linspace ( -1 . 33 , 1 . 32  ,  500  )  ; 
y^’  =  -.1083*xx  +  .2760; 

aa  =  linspace ( 1 . 32 , 1  -  245 , 500 )  ; 
bb  -  3.64’'aa  -  4.6718; 

cc  =  linspace ( -1 . 06 , 1 . 245 , 500 )  ; 
dd  =  .1533*cc  -.3344; 

ee  =  linspace ( -1 . 06 , -1 . 33 , 500 )  ; 
ff  =  -3.407*ee  -  4.1119; 

X  =  linspace (. 4099 ,. 57 14 , 500 )  ; 
y  =  -.14*x  +  .2785; 

a  =  linspace {, 5714 ,. 5273 , 500 )  ; 

L-  =  1.769»a  -  .8123; 

c  =  linspace( .5273,  .4102, 500)  ; 
d  =  -.2673*c  +  .26145; 

e  =  linspace (. 4102 ,  .4099,500); 

f  =  -232. 7*e  95.67 

split  =  .4098; 
xline  =  [split  split]; 
yline  =  [-1.32  1.32]  ; 

axis ( [-1.32  1.32  -.6  .5] )  ; 

hold  on;plot (xline, yline, 'g- ') ;  %  split 

plot  (xx, ^^-y ,  '  g .  '  ,  aa ,  bb,  '  g .  '  ,  cc,  dd,  '  g .  '  ,  ee,  f  f ,  '  g  .  '  )  ;  % 
quadrilateral 

plot (x, y , ' w. ' , a , b, ' w. ' , c , d, ' w. ' , e, f , ' w. ' ) ;  %  subquad 

t it le (' Subquad, Quadrilateral  and  Split  of.4098'); 

disp( 'about  to  print') 

print ; 

pause; clg; 

hold  off 

%while  s  0 

for  iter=l:3 

sprint f { 'Iter=%g'  ,iter) 


€ 


k 


9 


9 


9  i 
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for  n  =  1:500 

x0(n)  =  1.0  -  1.4*x(n)''2  +  y(n); 
yO (n)  =  . 3 'x (n) ; 
end 

%  1 

for  1  =  1:500 

a0{i)  =  1.0  -  1.4*a(i)''2  +  b{i); 
bO { i )  =  . 3  *a { i )  ; 
end 

%  2 

for  j  =  1:500 

c0(j)  =  1.0  -  1.4*c{j)''2  -h  d(j)  ; 
dO ( j )  =  . 3*c  (  j )  ; 


for  k  =  1:500 

eO (k)  =  1.0  -  1.4*e{k) "2  +  f  (k)  ; 
f0{k)  =  .3»e(k); 

end 

%  4 

axis{ [-1.32  1.32  -.6  .5] )  ; 

%  s  =  3  +  1  ; 

%subplot  (221 )  ; 

hold  on;plot (xline,yline, 'g- ' ) ;  %  split 

plot (xx, yy , ' g . ' , aa,  bb,  ' g . ' , cc, dd, ' g . ' , ee, f f , ' g . ' )  ;  % 
quadrilateral 

plot (xO , yO , ' w. ' , aO , bO , ' w. ' , cO , dO , ' w. ' , eO , f 0 , ' w. ' ) ;  % 
iterated  subquad 

title (' Iterated  Subquad, Quadrilateral  and  Split  of  .4098') 
pause ; 
hold  off 

%Prepare  for  next  iteration. 

x=x0 ;y=y0 ; a=a0 ; b=b0 ; c=c0 ; d=d0 ; e=e0 ; f =f 0 ; 

%  if  s==3 
print ; clg 
end 


APPIMDIX  B:  TYPICM.  LBFT  QtJASRILATBRAL  RUN 


The  following  is  a  typical  run  of  GRDC0MP6.M  which 
models  which  points  from  the  left  quadrilateral  give  the 
scune  binary  sequence  as  a  point  (W,2)  chosen  at  random  from 
the  left  quadrilateral.  The  point  (W,Z)  (+)  is  indicated 

with  a  long  arrow;  (u,v)  points  (x)  are  indicated  with 
shorter  arrows. 


Gnd  n«U  aad  Mac  Cham  M  Amctor 


Figure  55  Typical  run  left  quad  model:  before  iteration 


«) 
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Figure  56  Typical  run  left  quad  model:  Iterate  1 


Figure  58  Typical  run  left  quad  model:  Iterate  3 


•l  -W  0  OJ  I 

Figure  59  Typical  run  left  quad  model:  Iterate  4 


•  •  •  • 
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Flgur*  62  Typical  run  left  quad  model:  Iterate  7 


Flgur*  63  Typical  run  left  quad  model:  Iterate  8 


Flgurtt  64  Typical  run  left  quad  model:  Iterate  9 


Figure  65  Typical  run  left  quad  model:  Iterate  10 
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I  Ictracion 

Plgurtt  66  Grid  points  with  saune  sequence  after  n 
iterations 


Table  5  shows  the  count  of  points  for  iterates  1  through  10 
that  give  the  same  binary  sequence  as  (W,Z).  At  iterate  11 
(not  shown)  there  are  no  grid  points  that  give  the  same  binary 
sequence  as  (W, Z) . 


TABLE  5  NUMBER  OP  POINTS  THAT  GIVE  THE  SAME  BINARY  SEQUENCE 


NUMBER  OF  POINTS  THAT  GIVE  THE  SAME  BINARY  SEQUENCE 


ITERATE 

0 

■ 

2 

3 

■ 

5 

6 

■ 

8 

9 

10 

NO.  PTS 

314 

180 

150 

18 

3 

3 

3 

3 

3 

3 

3 
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Figure  67  (W,2)  and  chree  (u.v)  points  prior  to 
iteration 


Figure  68  (W,Z)  and  three  {u,v)'s:  Iterate  1 


Pigrurs  69  (W,Z)  and  three  (u,v)'s:  Iterate  2  (Notice 
that  the  (u,v)  points  are  essentially  superimposed). 


o(  Pt  Choieii  on  HA  lad  Mashiaf  S«qMc*  Pt  anal  Otv«t« 


Pigura  70  {W,Z)  and  three  {u,v)'s:  Iterate  3  (Notice 
that  the  (u,v)  points  are  essentially  superimposed) . 


Figure  71  (W,Z)  and  three  (u,v)'s;  Iterate  4  (Notice 
that  the  (u,v)  points  are  essentially  superimposed). 


Figure  72  (W,Z)  and  three  {u,v)'s:  Iterate  5  (Notice 
that  the  (u,v)  points  are  essentially  superimposed). 


VtovomM  of  Pi  OMMa  M  HA  lad  Mactac  Sa^yMM  Pt  MHil  Oivccf* 


.Vtovcmm  ft  CbOMi  os  HA  «ad  SaqutiK*  ft  hmU  Div«|« 


'J- 


OJ  1 


Figur*  74  (W,Z)  and  three  (u,v)'s:  Iterate  7 


Ma««Mm  o(  Pi  CboMB  oa  HA  Md  Si^onet  PI  Mdl  Oi*«fi 


Movcmm  of  Pt  CboM*  M  HA  wd  MMchuK  Saqwae*  Pt  aaal  Oi«ar|* 


APPBMDIZ  C:  6  TO  IS-TDPLX  DATA 

This  appendix  contains  the  counts  of  n-tuples  for  typical 
Henon  binary  sequences  for  n-tuples  of  length  6-16. 


Figure  78  Incidence  of  6-tuples  in  typical  Henon 
generated  binary  sequence 


Figure  79  Incidence  of  7 -tuples  in  typical  Henon 
generated  binary  sequence 
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Shir cad  Daclaal  Eaulvalanc 


rlgur*  80  Incidence  of  8-tuples  in  typical  Henon 
generated  binary  sequence 


Shlfced  DccImI  Equivalent 


rigure  81  Incidence  of  9 -tuples  in  typical  Henon 
generated  binary  sequence 


Shifted  OeclBiel  Equivalent 


Figure  82  Incidence  of  10 -tuples  in  typical  Henon 
generated  binary  sequence 


Figure  83  Incidence  of  11-tuples  in  typical  Henon 
generated  binary  sequence 


Figure  84  Incidence  of  12 -tuples  in  typical  Henon 
generated  binary  sequence 


Figure  85  Incidence  of  13 -tuples  in  typical  Henon 
generated  binary  sequence 
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Figurtt  86  Incidence  of  14 -tuples  in  typical  Henon 
generated  binary  sequence 


Figure  87  Incidence  of  15-tuples  in  typical  Henon 
generated  binary  sequence 


Figure  88  Incidence  of  16-tuples  in  typical  H6non 
generated  binary  sequence 
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APPIHDZX  D:  AN  ONRKALIZABLS  S-TOPLB 

As  we  saw  in  Chapter  IV,  :ust  because  a  particular  n-tuple 
is  not  found  in  a  typical  string  of  length  10‘  does  not  mean 
we  will  not  find  it  in  a  typical  string  of  length  10'  or 
longer.  The  only  way  to  show  that  an  n-tuple  is  unrealizable 
is  to  prove  it  analytically  as  we  have  done  for  the  sequence 
{1,1,  0,0}.  It  is  in  this  Appendix  that  we  prove  that  the 
sequence  (0,0, 0,0,0}  is  also  unrealizable.  This  specific  5- 
tupie  is  used  in  Chapter  V  in  order  to  support  a  conclusion  as 
to  why  certain  sequences  are  unrealizable. 

The  sequence  {0,0, 0,0,0}  is  one  of  the  four  5-tuples  that 
appeared  to  be  unrealizable  from  our  da^-a  in  addition  to  those 
5-tuples  that  included  the  sequence  (1,1, 0,0}.  These  5-tuples 
appear  in  Chapter  IV  Table  4.  Consider  the  5-tuple 
{0,0, 0,0,0^.  Suppose  this  binary  sequence  corresponds  to  the 
following  sequence  of  x  values:  {x.. ; ,  x..,  x..; ,  x.,.; ,  x..-, }  .  All 
these  X  values  must,  therefore,  be  at  most  .4098.  Referring 
to  equation  (1.3)  on  page  10  the  possible  solutions  for  x.  _ 
correspond  to  the  region  beneath  the  inverted  parabola  in 
Figure  89.  Similarly,  only  those  x,.-  and  x..  terms  less  than 
or  equal  to  .4098  are  considered  (refer  to  Figure  90) . 
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By  containing  the  trapping  region  in  a  rectangle  (see 
Figure  ?1)  we  see  how  the  possible  solution  values  can  be 
further  constrained. 


Nteii—ai  Btc— ■  PirampMi  Tit 


Figure  91  Trapping  region  in  minimum  area  rectangle 

The  minimum  and  maximum  x  values,  for  example,  are  -1.32  and 
1.33.  In  order  to  show  that  a  string  of  five  successive  zeros 
is  not  possible  we  consider  the  recurrence: 

x„.,  =  1  -  1.4x„.j2  +  .3x„.i  (D.l) 

where  x„  and  x^.i  are  implicit  in  the  equations.  To  optimize 
our  potential  of  finding  a  value  for  x„.j  at  most  x„e.d=-4098  we 
require  a  maximum  x„.2  value  and  a  minimum  x^.i  value.  Clearly, 


* 


W 


* 
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x„.;  can  be  no  less  than  -1.32.  However,  x„.2  is  not  limited 
only  by  .4098,  as  we  shall  see.  Consider  the  equation: 

x„.2  =  1  -  1.4x„,i^t-  .3x„. 

To  find  our  desired  maximum  x„.2  value,  we  must  use  a  minimum 
value  for  x^  and  a  maximum  value  for  x,..,  such  that 

x„.2  SI-  1.4(.4098)^  +  .3(-1.32). 


This  further  constrains  x„,2  to  a  value  no  greater  than 
approximately  .3689  (see  Figure  92).  (The  value  of  x.,.;  to 
eight  decimal  places  is  .36888954). 


sequence  (0,0, 0,0,0}  under  x„.i  and  x„.2  restrictions 
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Thus,  referring  to  equation  (D.l),  the  maximum  x,._  value  is  no 
more  than  .3689  and  the  minimum  x^,-.  value  is  -1.32  based  on 
the  confines  of  the  trapping  region.  It  follows  that 

X  >  1  -  1.4(  .3689)  .3(-1.32) 

x,..,  >  .41348871  (to  eight  decimal  places). 

The  value  of  x. .  .  corresponds  to  a  binary  value  of  1.  Thus 
{0,G, 0,0,0}  IS  not  realizable  under  this  symbolic  dynamics. 


4r 


142 


•  •  • 


J 


■*J 


« 


APPENDIX  S:  SUBORDINATE  PROGRAMS 

This  Appendix  includes  programs  which  are  necessary'  but 
secondary  to  the  programs  in  Appendix  A.  TV’d  of  the  program 
HENON  and  HENREAL)  are  essential  for  the  proper  operation  o 
the  APPENDIX  A  programs. 

function  [x,y]  a  henrttal(n,xO,yO) 

%  NOTE:  THIS  PROGRAM  IS  REQ'D  TO  RUN  ALL  GRDCOMP  PROGRAMS 
%  This  program  is  credited  in  its  entirety  to  the  author  of 
reference  13; it  has  been  used  by  permission, 
program  to  generate  n-length  real  sequences  based  or.  the  * 
Heno.n 

horseshoe  attractor.  initial  points  fit  into  the 
quadrilateral  of  convergence  described  in  [Kent 6]  . 
inputs : 

n  =  length  of  desired  sequence 
xO  =  initial  X  value 
yO  =  initial  y  value 
outputs : 

X  =  n  by  1  real  vector 
y  =  n  by  1  real  vector 

% 

X  =  zeros (n, 1 ) ; 
y  =  zeros (n, 1 ) ; 

X ( 1 )  =  xO ; 
y  ( 1 )  =  yO ; 

%routine  to  check  if  initial  points  are  valid 
A=[3.4074  1;-.1083  -1;  -3.64  1;  -.1562  1 ] ; 

B=[-4.1119  -.2760  -4.6718  -.3344]'; 
if  min ( ( A* [xO ;y0 ] )  >  B)  ==  0 

disp (' initial  point  outside  convergence  zone') 
return 

end 

%recursive  generation  of  points 
for  i  =  l:n-l; 

xli-t-l)  =  y(i)  +  1  -  1.4*x(i)"2; 
y  (  i  +  1 )  =  .3  *  X ( i ) ; 

end 

%plot (x,y, 'b. ' ) ; 


W 


Ar 
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%  PROBStTBQUAO.M 

%  PROGRAM  REQUIRES  HIOIREAL  PROGRAM 

%  This  program  finds  the  dynamic  probability  of  a  point 
falling  into  the 

%  exclusion  zone  compared  to  being  on  the  rest  of  the 
attractor . 

[x,y]  =  heareal ( 100000 , 0 , 0 ) ; 

a  =  [232.7  1;  -.14  -1;  -1.769  1;  .2673  1] ; 

b  =  [95.6  -.29  -.8  .25] ' ; 

-%z  =  zeros { length (x)); 


f  cr  1  = 
subquad 


1 : length ( x) ; 

if  min { (a*[x{i) ;y(i) ] ) 

z{i)  =  1; 
else  z ( i )  =  0 ; 
end 


end 

t  =  f ind ( z  =  =  0  )  ; 
tt  =  length ( t )  ; 


>  b) 


0  %  if  outside 


prcbabilityllO  =  tt/length(x) 
pause; 

subplot (2,1,1) ,plot{x(t) ,y(t) , 'w.') 
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function  [xaat,  ynat]  «  grdconpS ( sp ) 

ONLY  DIFF  BTWN  4&5  IS  THAT  4  CHOOSES  W,  Z  FROM  THE  LEFT 
UAD  THEN  ENSURES 

THIS  POINT  IS  ON  THE  ATTRACTOR,  5  CHOOSES  FROM  THE  RIGHT 
UAD  AT  RANDOM 

THEN  ENSURES  THAT  THIS  POINT  IS  ON  THE  ATTRACTOR. 

This  function  takes  a  point (W, Z) from  the  quadr  and 
terates  the  point  n 

times  in  order  to  ensure  the  point  is  on  the  Henon 
1 1  ractor ( ensuring 

that  the  nth  iterate  is  in  the  left  quadrant) .  It 
nitiates  a  grid  field 

based  on  a  certain  spacing,  than  iterates  those  grid  field 
cints  that 

match  the  binary  string  of  (wnew, znew) .  As  long  as 
the  strings  match,  points  are  iterated  using  the  Henon 
ecurrence.  Only  that 

those)  points  that  completely  match  are  finally  plotted 
ccording  to  the 
f cl lowing  scheme; 

g.  original  grid  field  points 

w*  original  (W,Z) 

gx  original  xvec,yvec  point (s)  which  matches  binary 

tring  of  W, Z 

crmat  long 

We  also  show  how  W,Z  and  xvec ( indices ), yvec ( indices ) 
walk"  to  the  Henon 

attractor  at  each  iteration  (with  the  attractor  on 
creen)  . 

=linspace (-1.33,1.32,500); 

=-.1083*0  +  .276; 

=  linspace (1.32, 1.245, 500)  ; 

=3.64»u  -4.6718; 

=linspace (-1.06,1.245,500) ; 

=.1533*g  -.3344; 

=linspace (-1.06,-1.33,500); 

=-3.407*e  -  4.1119; 
split  =  .4098; 

fplx  =  .6314;  fply  =  .1894; 

fp2x  =  -1.1314;  fp2y  =  -.3393; 
a  =  [  -.1083  -l;-3.64  1;  -.1562  1]; 
b  =  [  -.2760;  -4.6718;  -.3344]; 

C  =  [ split : sp : 1  .  32 ]  ; 
d  =  [ - . 6 : sp :  . 5 ]  ; 

[x,  y]  =  henreal ( 500 , - 1 . 0 , - . 25 )  ; 

xline  =  [split  split];  yline  =  [-1.32  1.32]; 

lenc  =  length(c);  lend  =  iength(d); 

xvec  =  zeros ( 1 , lenc*lend) ;  yvec  =  zeros (1, lenc*lend/  ; 

k=G; 

for  p  =  1 : lenc 
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for  r  =  l:lend 
k  =  k  +  1  ; 

%  This  next  section  puts  everything  outside  of  the  quad  to 
;  -10, -10) 

if  max{(a  *  [c{p);d(r)])  <  b)  ==  0 
xvec{k)  =  c(p);  yvec (k) =d{r) ; 
else 

xvec(k)  =  -10;  yvec(k)  =  -10; 

end 

end 

end 

xvec  =  reshape (xvec, lenc*lend, 1) ; 
yvec  =  reshape (yvec , lenc*lend, 1 ) ; 

m  =  findixvec  ~=  -10  I  yvec  ~=  -10); 
xvec  =  xvec (m) ; 
yvec  =  yvec (m) ; 

num  =  rand*length (xvec ) ; 
index  =  fix (num)  +  1; 

w  =  xvec(index); 
z  =  yvec ( index) ; 

wattr  =  w;  %  This  preserves  the  values  of  W,Z 

zattr  =  z; 

for  j  =  1:40  % 

ensures  wattr, zattr  is 

wattr(j+l)  =  1.0  -  1.4*wattr ( j ) ^2  +  zattr(j);%  on  the 
attractor 

zattr(j+l)  =  .3*wattr(j); 

if  j  >  20  Sc  wattr(j+i)  >  split,  break; 
end 

end 

wnew  =  wattr(j+l); 
znew  =  zattr (j+1) ; 
wattr  =  wnew; 
zattr  =  znew; 

XX  =  xvec;  %  preserves  the  values  of 

xvec  x  yvec 
yy  =  yvec; 

axis( [-1.32  1.32  -.6  .5]  )  ; 

indices  =  1 : length (xx) ;  %  (1  2  3  4  .  .  .] 

1  =  1; 

qvec  =  length(xx); 

q  =  indices; 
lenq  =  length (xx) 
hold  on; 

plot (xvec, yvec, 'g. ' , wattr, zattr, 'w+' , xline,yline, 'w' ) ; 

plot (o,s, 'w.',u,v, 'w.',e,f,'w.',g,h,'w.', fplx, fply, 'o' , fp2x, 
fp2y, 'o' ) ; 
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hold  off; 

title('Grid  Field  and  Point  Chosen  on  Attractor'); 

print ; 
pause; clg • 

while  length (q)  >  1; 

axis( [-1.32  1.32  -.6  .5] )  ; 

xO  =  yy  +  ones (size (yyO  )  -  1.4*xx."2; 

yO  =  .3*xx; 

wO  =  znew  +  1  -  1.4*wnew."2; 
zO  =  .3*wnew; 

if  wO  >=  split 

qtemp  =  find(xO  >  split  I  xO  ==  split); 
else 

qtemp  =  findixO  <  split); 
end; 

if  length ( qtemp)  ==  0, break; 

else  q  =  qtemp; 

end; 

i  =  i  +  1;  %  counts  iterates  where  at  least 

one  point  matches 

qvec  =  [qvec , length ( q) ]; %  plot  of  how  #  with  same  binary 
sequence 

lenq  =  length (q) 

indices  =  indices (q);  %  decreases  each  time  through  loop 


XX  =  x0(q);  %  preserves  to  next  iterate  those 

matching  values  of 

yy  =  y0(q);  %  xO  and  yO 

wnew  =  wO ; 
znew  =  zO; 
pause ; 
hold  on; 

plot (XX, yy,  ' g . ' , wO , zO , ' w+ ' , xline,yline,  ' w  '  )  ; 


plot(o,s,  'w.  ',u,v,  'w.',e,f,  'w.',g,h,  'w.',fplx,fply,  'o',fp2.x, 
fp2y,  'o' )  ; 

title ( 'Movement  of  Points  in  Grid  with  Same  Binary 
Sequence ' ) ; 
print ; 

pause;  clg;  hold  off; 


end 

axi3{[l  2  3  4]);  axis; 

vec  =  ( 1/2 )."(( 1 : length (qvec) ) -1 ) ; 

plot ( ( 1 : length (qvec ) ) -1 , qvec, ' r ' , ( 1 : length (qvec) ) -1 , qvec ( 1 ) * 
vec , ' b ' ) ; 

title ( 'Decrease  in  #  of  Grid  Points  with  Same  Binary 
Sequence  vs  l/2^n') 
print ; 
pause ; clg ; 


•  • 


•  • 
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%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  original 
xvec  and  yvec . 

axis( [-1.32  1.32  -.6  .5]  )  ; 
hold  on; 

plot(xline,yline,'w',x,y,'b.',o,s,'w.',u,v,'w.',e,f,'w.',g,h 
, 'w. ' ) ;hold  on; 

plot  (xvec  {  indices )  , yvec  ( indices )  ,  '  gx'  ,  wattr ,  zattr ,  '  w-*- '  )  ; 
plot ( fplx, f ply , ' o ' , fp2x, fp2y , 'o'); 
title ( 'wattr, zattr  &  Grid  Pt  with  Same  Binary  Sequence'); 
print ; 

pause ;clg; hold  off; 
xnew  =  xvec ( indices ) ; 

ynew  =  yvec ( indices ) ;  %  now  we  show  a  plot  of  how 
wattr, zattr  and  the  point 

%  that  generates  the  same  binary  sequence 

xnew, ynew 

%  walk  around  the  attractor 

xnewt  =  xnew; 
ynewt  =  ynew; 
wattrt  =  wattr; 
zattrt  =  zattr; 
for  n  =  1 ; 1  -  1 

xOt  =  1.0  -  1.4*xnewt.^2  +  ynewt; 
yOt  =  .3*xnewt; 

wOt  =  1.0  -  1 . 4*wattrt . ''2  +  zattrt; 
zOt  =  . 3 'wattrt ; 
wattrt  =  wOt; 
zattrt  =  zOt; 
xnewt  =  xOt; 
ynewt  =  yOt; 

axis( [-1.32  1.32  -.6  .5] )  ; 
hold  on; 

plot(x,y,'b.', xline, yline, 'w-',o,s, 'w.',u,v, 'w.',e,f, 'w.',g, 
h ,  '  w .  '  )  ; 

plot (fplx, fply , 'o' , fp2x, fp2y , 'o'); 
plot (xOt , yOt , ' gx ' , wOt , zOt , ' w+ ' ) ; 
title ( 'Movement  of  Pt  Chosen  on  HA  and  Matching  Sequence  Pt 
until  Diverge'); 

print ; 
pause; clg; 

%plot (xOt ,y0t , ' ix' ,w0t, zOt, ' i+' ) ; 
end 

hold  o f f ; 
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%  function  [xnat,  yaat]  ■  grdcoac>6  (sp) 

%  ONLY  DIFF  BTWN  4&6  IS  THAT  4  USED  WATTR, ZATTR  5  CHOOSES  AT 
RANDOM  FROM 

%  THE  LEFT  QUADR  THEN  FINDS  WATTR, ZATTR. 

%  This  function  takes  a  point (W, Z) from  the  quadr  and 
Iterates  the  point  n 

%  times  in  order  to  ensure  the  point  is  on  the  Henon 
attractor ( ensuring 

%  that  the  nth  iterate  is  in  the  left  quadrant) .  It 
initiates  a  grid  field 

%  based  on  a  certain  spacing,  then  iterates  those  grid  field 
points  that 

%  match  the  binary  string  of  (wnew, znew) .  As  long  as 
%  the  strings  match,  points  are  iterated  using  the  Henon 
recurrence.  Only  that 

%  (those)  points  that  completely  match  are  finally  plotted 
according  to  the 
%  following  scheme: 

%  g.  original  grid  field  points 

%  w+  original  (W,Z) 

%  gx  original  xvec,yvec  point (s)  which  matches  binary 

string  of  W,Z 

format  long 

%  We  also  show  how  W, Z  and  xvec ( indices ), yvec ( indices ) 

"walk"  to  the  Henon 

%  attractor  at  each  iteration  (with  the  attractor  on 
screen) . 

o  =  linspace (-1.33,1.32,500)  ; 
s=-. 1083*0  +  ;276; 
u=linspace (1.32,1.245,500)  ; 
v=3.64*u  -4.6718; 
g=linspace(-1.06, 1.245, 500) ; 
h=.1533*g  -.3344; 
e=l inspace (-1.06,-1.33,500)  ; 
f=-3 .407*6  -  4.1119; 
split  =  .4098; 

fplx  =  .6314;  fply  =  .1894; 

fp2x  =  -1.1314;  fp2y  =  -.3393; 
a  =  [3.4074  1;  -.1083  -1;  -.1562  1]; 
b  =  [-4.1119;  -.2760;  -.3344]  ; 

c  =  [ -1 . 32 : sp : split ] ; 
d  =  [ - . 6 : sp :  .  5 ]  ; 

[x,  y]  =  henreal ( 500 , -1 . 0 , -  .  25 )  ; 

xline  =  [split  split];  yline  =  [-1.32  1.32]; 

lenc  =  length(c);  lend  =  length(d); 

xvec  =  zeros ( 1 , lenc*lend) ;  yvec  =  zeros ( 1 , lenc* lend) ; 

k=0; 

for  p  =  l:lenc 

for  r  =  l:lend 
k  =  k  +  1  ; 
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%  This  next  section  puts  everything  outside  of  the  quad  to 
(-10,-10) 

if  max((a  *  [c(p);d(r)])  <  b)  ==  0 
xvec(k)  =  c(p);  yvec (k) =d ( r ) ; 
else 

xveclk)  =  -10;  yvec(k)  =  -10; 

end 

end 

end 

xvec  =  reshape (xvec, lenc’lend, 1 > ; 
yvec  =  reshape (yvec , lenc*lend, 1 ) ; 

m  =  find(xvec  ~=  -10  I  yvec  -=  -10); 
xvec  =  xvec (m) ; 
yvec  =  yvec (m) ; 

num  =  rand*length (xvec) ; 
index  =  fix (num)  +  1; 

w  =  xvec (index); 
z  =  yvec ( index) ; 

wattr  =  w;  %  This  preserves  the  values  of  W, Z 

zattr  =  z; 

for  j  =  1:40  % 

ensures  wattr, zattr  is 

wattr(j+l)  =  1.0  -  1 .4*wattr  ( j ) ''Z  +  zattr(j);%  on  the 
attractor 

zattr(j+l)  =  .3*wattr(j); 

if  j  >  20  &  wattr(j+l)  <  split,  break; 
end 

end 

wnew  =  wattr(i+l); 
znew  =  zattr ( j  +1 ) ; 
wattr  =  wnew; 
zattr  =  znew; 

XX  =  xvec;  %  preserves  the  values  of 

xvec  &  yvec 
yy  =  yvec; 

axis( [-1.32  1.32  -.6  .5]  )  ; 

indices  =  1 : length (xx) ;  %  [1  2  3  4  .  .  .] 

i  =  1; 

qvec  =  length(xx); 

q  =  indices; 
lenq  =  length (xx) 
hold  on; 

plot (xvec , yvec , ' g . ' , wattr, zattr, ' w+ ' , xline, yline, ' w' ) ; 

plot(p,s, 'w.',u,v, 'w.',e,f, 'w.',g,h, 'w.',fplx,fply, 'o',fp2x, 
fp2y, 'o' ) ; 
hold  off; 

title('Grid  Field  and  Point  Chosen  on  Attractor'); 


print ; 
pause ; clg ; 

while  length (q)  >  1 ; 

axis( [-1.32  1.32  -.6  .5] )  ; 

xO  =  yy  +  ones  (size  (yy)  )  -  1.4’xx.''2; 

yO  =  .3*xx; 

wO  =  znew  +  1  -  1 . 4*wnew. "2 ; 
zO  =  .3*wnew; 

if  wO  >=  split 

qtemp  =  find(xO  >  split  I  xO  ==  split); 
else 

qtemp  =  findtxO  <  split); 
end; 

if  length (qtemp)  ==  0, break ; 

else  q  =  qtemp; 

end; 

i  =  i  +  1;  %  counts  iterates  where  at  least 

one  point  matches 

qvec  =  [qvec ,  length  ( q)  ];  %  plot  of  how  #  with  same  binar^^ 
sequence 

lenq  =  length(q) 

indices  =  indices (q);  %  decreases  each  time  through  loop 

XX  =  xO(q);  %  preserves  to  next  iterate  those 

matching  values  of 

yy  =  yO (q) ;  %  xO  and  yO 

wnew  =  wO ; 
znew  =  zO; 
pause; 
hold  on; 

plot ( XX, yy , ' g . ' , wO , zO , ' w+ ' , xline,yline, ' w' ) ; 

plot(o,s, 'w.',u,v, 'w.',e,f, 'w.',g,h, 'w.',fplx,fply, 'o',fp2x, 
fp2y, 'o'); 

title ( 'Movement  of  Points  in  Grid  with  Same  Binary 
Sequence ' ) ; 
print ; 

pause;  clg;  hold  off; 


M> 


Ar 


•  • 


end 

axi3([l  2  3  4]);  axis; 

vec  =  ( 1/2  ).''(( 1 :  length  (qvec)  ) -1 )  ; 

plot ( ( 1 ; length (qvec ) ) -1 , qvec , ' r ' , ( 1 : length (qvec ) ) -1 , qvec ( 1 ) * 
vec ,'b'); 

title (' Decrease  in  #  of  Grid  Points  with  Same  Binary 
Sequence  vs  l/2''n') 
print 
pause ; clg; 

%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  grid  point  the 


_ _ . 


M. 
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%  surviving  iterate  corresponds  with  respect  to  the  original 
xvec  and  yvec . 

axis( [-1.32  1.32  -.6  .5] ) ; 
hold  on; 


plot (xline,yline, 'w',x,y, 'b. ',o,s, 'w. ',u,v, 'w. ',e,f, 'w. ' ,g,h 
,  '  w .  '  )  ; 

plot  (xvec  (  indices)  ,yvec  (indices)  ,  'gx'  ,  wattr,  zattr,  'w-r '  )  ; 
plot (fplx, fply, 'o' , fp2x, fp2y , 'o'); 
title ( 'wattr , zattr  &  Grid  Pt  with  Same  Binary  Sequence'); 

print ; 
pause; clg; 
hold  o f f ; 

xnew  =  xvec ( indices ) ; 

ynew  =  yvec ( indices ) ;  %  now  we  show  a  plot  of  how 
wattr, zattr  and  the  point 

%  that  generates  the  same  binary  sequence 

xnew, ynew 

%  walk  around  the  attractor 

xnewt  =  xnew; 
ynewt  =  ynew; 
wattrt  =  wattr; 
zattrt  =  zattr; 
for  n  =  l:i  -  1 

xOt  =  1.0  -  1 . 4*xnewt . ''2  +  ynewt; 
yOt  =  .3*xnewt; 

wOt  =  1.0  -  1 . 4 *wattrt . ^2  +  zattrt; 
zOt  =  .3*wattrt; 
wattrt  =  wOt; 
zattrt  =  zOt; 
xnewt  =  xOt; 
ynewt  =  yOt; 

axis! [-1.32  1.32  -.6  .5] ) ; 
hold  on; 


•  • 


plot (x, y , ' b . ' , xline,yline, 'w-',o,s, 'w.',u,v, 'w.',e,f, 'w.',g, 
h ,  '  w .  '  )  ; 

plot (fplx, fply, 'o' , fp2x, fp2y , 'o') ; 
plot (xOt ,y0t , ' gx ' , wOt , zOt , ' w+ ' ) ; 
title { 'Movement  of  Pt  Chosen  on  HA  and  Matching  Sequence  Pt 
until  Diverge'); 

print ; 
pause; clg; 

plot(x0t,y0t, 'ix',w0t,z0t, 'i+'); 
end 

hold  off; 


%  function  [xmat,  ymat]  «  proof 5 (apcrs, incr, apfn) 

%  HERE  WE  USE  SPCRS  TO  CHOOSE  A  POINT  AT  RANDOM  THEN  ENSURE 
IT  IS  ON  THE 

%  ATTRACTOR.  IT  IS  THEN  USED  AS  THE  COMPARISON  POINT  FOR  ALL 
FURTHER  SP  VALUES 

%  ONLY  DIFF  BTWN  5i6  IS  THAT  6  CHOOSES  W, Z  FROM  THE  LEFT 
QUADR  THEN  ENSURES 

%  THIS  P-^INT  IS  ON  THE  ATTRACTOR,  5  CHOOSES  FROM  THE  RIGHT 
QUADR . 

%  A  NEW  WATTR.ZATTR  ARE  FOUND  EACH  TIME  SP  CHANGES. 


%  This  function  takes  a  point :w, Z ; from  the  quadr  and 
Iterates  the  point  n 

%  times  in  order  to  ensure  the  point  is  on  the  Henon 
attractor! ensuring 

%  that  the  nth  iterate  is  in  the  left  quadrant) .  It 
initiates  a  grid  field 

%  based  on  a  certain  spacing,  then  iterates  those  grid  field 
points  that 

%  match  the  binary  string  of  (wnew, znew) .  As  long  as 
%  the  strings  match,  points  are  iterated  using  the  Henon 
recurrence.  Only  that 

%  (those)  points  that  completely  match  are  finally 
considered . 
split  =  .4098; 
sp  =  spcrs; 

a  =  [  -.1083  -l;-3.64  1;  -.1562  1  ]  ; 
b  =  [  - .2^60;  -4.6718;  -.3344] ; 
c  =  [ spl It : sp : 1  .  32 1  ; 
d  =  [ - . 6 : sp :  .  5 ]  ; 

lenc  =  length(c);  lend  =  length(d); 

xvec  =  zeros ( 1 , lend* lenc) ;  yvec  =  zeros ( 1 , lend* lenc ) ; 

k=0; 

for  p  =  1 : lenc 

for  r  =  1 : lend 

k=k+l ; 

1  f  max  ;;a  *  [c(p);d!r)])  <  b)  0 
xvec(k)  =  c(p);  yvec ( k ) =d ( r ) ; 
else 


end 

end 

xvec  = 
yvec  = 
m  = 
xvec  = 
yv0c  = 


end 


xvec(k)  =  -10; 


reshape (xvec, lenc* lend, 1) ; 
reshape (yvec , lenc* lend, 1 ) ; 
find{xvec~=  -10  I  yvec--=  -10); 
xvec (m) ; 
yvec (m) ; 


yvec (k) 


-10; 


num 

index 


rand*length (xver )  ; 
fix (num)  +  1; 


%  This  preserves  the  values  of 


w  =  xvec ( index) ; 
z  =  yvec (index); 
wattr  =  w; 

W,  Z 

zattr  =  z; 
for  j  =1:40 

%ensures  wattr, zattr  is 

wattr(j+l)  =  1.0  -  1 . 4*wattr ( j ) "2  +  zattr(j);%  on  the 
attractor 

zattr (j+1)  =  .3*wattr(j); 

if  D  >  20  Sc  wattr{3  +  l)  >  split,  brealc; 
end 

end 

wnew  =  wattr(j+l); 
znew  =  zattr ( j  +  1 )  ; 

for  sp  =  spfn : incr : spcrs ,  %  MAJOR  OUTER  LOOP 

sp  =  sp 

a  =  [  -.1083  -l;-3.64  1;  -.1562  1]; 
b  =  [  -.2760;  -4.6718;  -.3344]; 
c  =  [ split ; sp : 1 . 32 ]  ; 
d  =  [ - . 6 : sp : . 5 ] ; 

lenc  =  length(c);  lend  =  length(d) ; 

xvec  =  zeros ( 1 , lend* lenc ) ;  yvec  =  zeros ( 1 , lend’lenc ) ; 

k=0; 

for  p  =  l:lenc 

for  r  =  l:lend 

k  =  )<+l  ; 

i f  max ((a  *  [c{p);d(r)])  <  b)  ==  0 
xvec(k)  =  c(p);  yvec (k) =d{r) ; 
else 

xvec(k)  =  -10;  yvec{k)  -10; 

end 

end 

xvec  =  reshape (xvec, lenc*lend, 1) ; 
yvec  =  reshape (yvec , lenc*lend, 1 ) ; 

m  =  find(xvec~=  -10  I  yvec-=  -10); 
xvec  =  xvec(m); 
yvec  =  yvec (m) ; 

wattr  =  wnew; 
zattr  =  znew; 

XX  =  xvec;  %  preserves  the  values  of 

xvec  Sc  yvec 
yy  =  yvec; 

indices  =  1 : length (xx) ;  %  [1234  .  .  .] 

i  =  1; 

qvec  =  length (xx); 
q  =  indices; 
lenqvec  =  length(xx); 
while  length (q)  >  1; 


xO  =  yy  +  ones  ( size  (yy)  )  -  1.4*xx.''2; 
yO  =  .3*xx; 

wO  =  znew  +  1  -  1 . 4*wnew . ^2 ; 
zO  =  .3*wnew; 

if  wO  >=  split 

qtemp  =  find(xO  >  split  I  xO  ==  split); 
else 

qtemp  =  find(xO  <  split); 
end; 

if  length ( qtemp )  ==  0, break; 

else  q  =  qtemp; 

end; 

1=1+  l;%counts  iterates  where  >  one  point  matches 
qvec  =  [qvec,  length(q)  ]  ;  %  plot  of  how  #  with  same  binary' 
sequence 

lenq  =  length(q); 

indices  =  indices (q);  %  decreases  each  time  through  loop 

XX  =  x0(q);  %  preserves  to  next  iterate  those 

miatching  values  of 

yy  =  y0(q) ;  %  xO  and  yO 

wnew  =  wO ; 
znew  =  zO; 

end 

%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  original 
xvec  and  yvec . 

xnew  =  xvec ( indices ) ; 
ynew  =  yvec ( indices )  ; 
xnewt  =  xnew; 
ynewt  =  ynew; 
wattrt  =  wattr; 
zattrt  =  zattr; 
for  n  =  1 : i  -  1 

xOt  =  1.0  -  1 . 4*xnewt  . ''2  +  ynewt; 
yOt  =  .3*xnewt; 

wOt  =  1.0  -  1 . 4*wattrt . ''2  +  zattrt; 
zOt  =  .3*wattrt; 
wattrt  =  wOt; 
zattrt  =  zOt; 
xnewt  =  xOt; 
ynewt  =  yOt; 
end 

nvec  =  [nvec,n]; 
spvec  =  [spvec,sp]; 
clear  xvec;  clear  yvec; 

end  %  END  TO  MAJOR  OUTER  LOO? 

axis  ; 

plot ( spvec , nvec , ' w* ' ) 

title('SP  VS.  ITERATIONS  UNTIL  DIVERGENCE'); 
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%  fiinction  [xaat,  yaat]  •  proof  6  (spcrs,  Incr,  spfn) 

i  HERE  WE  USE  SPURS  TO  FIND  A  POINT  THEN  ENSURE  IT  IS  ON  THE 
ATTRACTOR . 

%  THIS  POINT  IS  THEN  USED  AS  A  COMPAP.ISON  POINT  FOR  MANY  SP 
VALUES . 

%  ONLY  DIFF  BTWN  4i6  IS  THAT  4  USED  WATTR, ZATTR  6  CHOOSES  AT 
RANDOM  FROM 

%  THE  LEFT  QUADR  THEN  FINDS  WATTR, ZATTR. 

%  This  function  takes  a  point  .W, z ; from  the  quadr  and 
Iterates  the  point  n 

%  times  in  order  to  ensure  the  point  is  on  the  Henon 
at  t  ract or ' ensuring 

i  that  the  nth  iterate  is  in  the  left  quadrant) .  It 
initiates  a  grid  field 

%  based  on  a  certain  spacing,  then  iterates  those  grid  field 
points  that 

i  match  the  binary  string  of  {wnew, znewi .  As  long  as 
%  the  strings  match,  points  are  iterated  using  the  Henon 
recurrence.  Only  that 

%  : those'  points  that  completely  match  are  finally 
considered . 


sp- 


=  .4098; 


-.1083  -1  ;  -.1562  1]  ; 
-.2760;  -.3344]; 

split  ]  ; 


sp  =  spcrs; 
a  =  [3 .4074  1 
b  =  [-4.1119; 
c  =  [ -1 . 32  :  sp 
d  =  [ - . 6  :  sp ;  .  5 ]  ; 

lenc  =  length  c);  lend  =  length(d); 

xvec  =  zeros ( 1 , lend* lenc ) ;  yvec  =  zeros ( 1 , lenc* lend ) ; 

k  =  0  ; 

for  p  =  1 : lenc 

for  r  =  1 : lend 
k=k^l; 

if  max ((a  *  [c(p);d{r)])  <  b)  ==  0 
xvec(k)  =  c'p);  yvec ( k; =d ( r ) ; 
else 

xvec(k)  =  -10;  yvec(k)  -10; 

end 

end 

end 

.xvec  =  reshape  (xvec ,  lenc* lend,  1 )  ; 
yvec  =  reshape (yvec , lenc* lend, 1 ) ; 

m.  =  find(xvec~=  -10  I  yvec~=  -10); 
xvec  =  xvec(m); 
y'/ec  =  yvec  (m)  ; 


•  • 


num  =  rand*length(xvec) ; 
.ndex  =  fix  mum)  +  1  ; 


w  =  xvec ( index); 
z  =  yvec (index); 
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%  This  preserves  the  values  of  W,  Z 


water  =  w; 
zattr  =  z; 
for  ]  =  1:40 
%  ensures  ..attr,  zattr  is 

wattri]+l)  =  1.0  -  1 . 4 ’water  ( :  ) ''2  zattri^);!  on  the 
attractor 

zattr  1 1  +  1)  =  .3’wattr!j}; 

if  ]  >  20  i  water  !0:  +  l)  <  split,  break; 
end 

end 

-  wnew  =  wactr(i+l); 
znew  =  zattr ( 3  +  1 j  ; 

for  sp  -  spf n : incr : spers ,  %  MAJOR  OUTER  LO' 

sp  =  sp 

a  =  [3.4074  1;  -.1083  -1;  -.1562  1]; 
b  =  [-4.1119;  -  2760;  -.3344]  ; 

c  =  [ -1 . 32 : sp : split ]  ; 
d  =  [ - . 6 : sp :  .  5 ]  ; 

lenc  =  length(c);  lend  =  length(d); 

xvec  =  zeros ( 1 , lend* lenc ) ;  yvec  =  zeros ( 1 , lenc’ lend ) ; 

k30; 

for  p  =  1 : lenc 

for  r  =  l:lend 
k=k+l ; 

if  max ((a  ’  [ c (p ) ; d ( r ) ] )  <  b)  ==  0 
xvec(k)  =  c(p);  yvec ( k) =d ( r ) ; 
else 

xvec(k)  =  -10;  yvec(k)  -10; 

end 

end 

end 

xvec  =  reshape (xvec, lenc’lend, 1) ; 
yvec  =  reshape (yvec, lenc’lend, 1) ; 

m  =  find(xvec-=  -10  I  yvec-'=  -10); 
xvec  =  xvec(m); 
yvec  =  yvec (m) ; 
wattr  =  wnew; 
zattr  =  znew; 

XX  =  xvec;  %  preserves  the 

values  of  xvec  &  yvec 
yy  =  yvec; 

indices  =  1 : length (xx) ;  %  [1234  .  .  .] 

1  =  1; 

qvec  =  length(xx); 
q  =  indices; 
lenqvec  =  length(xx); 
while  length (q)  >  1; 

xO  =  yy  +  ones  ( size  (yy)  )  -  1.4*xx.''2; 
yO  =  .3*xx; 

wO  =  znew  +  1  -  1 . 4*wnew. ''2  ; 
zO  =  . 3 ’wnew; 


•  •  • 


•  .• 
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if  wO  >=  split 

qtemp  =  findixO  >  split  I  xQ  ==  split); 
else 

qtemp  =  find(xO  <  split); 
end; 

if  length  (qtemp)  ==  0,  break.  ; 
else  q  =  qtemp; 
end  ; 

1  =  1  1 ;  %  counts  iterates  where  at 

least  one  point  matches 

qvec  =  [ qvec , length ( q) ] ;  %  plot  of  how  #  with  same  binarA' 
sequence 

lenq  =  length (q) ; 

indices  =  indices (q);  %  decreases  each  time  through  loop 

XX  =  x0(q);  %  preserves  to  next  iterate  these 

matching  values  of 

'/Y  -  yC'q)  ;  %  xO  and  yO 

wnew  =  wO ; 
znew  =  zO; 

end 

%  In  this  way,  we  need  the  index  key  'indices'  to  tell  to 
which  grid  point  the 

%  surviving  iterate  corresponds  with  respect  to  the  original 
xvec  and  yvec . 

xnew  =  xvec ( indices ) ; 
ynew  =  yvec ( indices ) ; 
xnewt  =  xnew; 
ynewt  =  ynew; 
wattrt  =  wattr; 
zattrt  =  zattr; 
forn=l:i-l 

xOt  =1.0  -  1.4*xnswt.^2  +  ynewt; 
yOt  =  . 3 ’xnewt; 

wOt  =  1.0  -  1 . 4*wattrt . ''2  +  zattrt; 
zOt  =  . 3 ’wattrt; 
wattrt  =  wOt; 
zattrt  =  zQt; 
xnewt  =  xOt; 
ynewt  =  yOt; 
end 

.nvec  =  [nvec,n]; 
spvec  =  [spvec,sp]; 
clear  xvec; clear  yvec; 

end  %  END  TO  MAJOR  LOOP 

axis  ; 

plot ( spvec , nvec , ' w* ' ) 

title ('SP  VS.  ITERATIONS  UNIIL  DIVERGENCE'); 


%  END  TO  MAJOR  LOOP 


cx>  (D  o\“  ffi  ov'  X  (T)  lu  ;i'  "'f’  w  x  ■->*(' 


function  x  «  h«non(n,xO,yO) 

%  This  program  is  credited  to  ._ne  author  of  reference  13;  i 
%  has  been  used  by  permission. 

%  function  X  =  henon ' n , xO , yO ) 

%program  to  generate  n-length  binary  sequences  based  on  the 
% Henon 

%hcr3eshoee  attractor.  i.nitiai  points  are  checked  %against 

quadrilateral  of  convergence.  (Kenon  1976) 
inputs  : 

n  =  length  of  desired  sequence 
xO  =  initial  x  value 
yO  =  initial  y  value 

output s  : 

X  =  n  by  1  binaiT’  vector 

'1  =  xO ; 

:  1 '  =  yO ; 

plit  =  .4098;  %median  x-value  of  henon  attractor 

routine  to  check  if  initial  points  are  valid 
=[3.4074  1;-.1083  -1;  -3.64  1;  -.1562  1]; 

=[-4.1119  -.2760  -4.6718  -.3344]'; 
f  min  '  [xO  ;y0  n  >  B)  ==  0 

disp (' initial  point  outside  convergence  zone') 
return 
nd 

;2:n)  =  zeros' (n- 1 ,  1 )  ;  %vectors  are  preallocated  here  to 

%save  time 

'.2;n;  =  zeros  (n-1 ,  1 )  ;  %in  case  initial  point  is  outside 

%of  zone 

recursive  generation  of  points 
or  1  =  l;n; 

x(i  +  l)  =  y(i)  +  1  -  1.4*x(i)''2; 
y(i+l)  =  .3  *  x(i); 

convert  previous  point  to  binary 


if  x(i) 

<=  split 

X  {  i  )  =  0  ; 

else 

x(i)  =  1; 

end 

nd 

miinor  housekeeping  to  dump  the  last  term 

X  =  X  ( 1 ;  n )  ; 
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